From 318a254e9162677dd0e634556abed0ec39508ad9 Mon Sep 17 00:00:00 2001 From: Michael Zingale Date: Fri, 1 Nov 2024 13:49:43 -0400 Subject: [PATCH 1/4] update to 24.11 --- external/Microphysics | 2 +- external/amrex | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/external/Microphysics b/external/Microphysics index 533f712d2d..1753bb908e 160000 --- a/external/Microphysics +++ b/external/Microphysics @@ -1 +1 @@ -Subproject commit 533f712d2d217ea4f4063f3a809eb16a839c9ef6 +Subproject commit 1753bb908e2a966e66a3a79872f4c496c1df9fd3 diff --git a/external/amrex b/external/amrex index aba5f01be1..4a6e7c87e1 160000 --- a/external/amrex +++ b/external/amrex @@ -1 +1 @@ -Subproject commit aba5f01be1f0dbd93600671f6d26df228b7cca12 +Subproject commit 4a6e7c87e1494d10893b36f0d7ce6bf9df19fe0d From 32e14bff48d1bc6a11f52b422bcffaa47e8eb5b6 Mon Sep 17 00:00:00 2001 From: Michael Zingale Date: Fri, 1 Nov 2024 13:50:26 -0400 Subject: [PATCH 2/4] update changes for 24.11 (#2985) --- CHANGES.md | 21 +++++++++++++++++++++ 1 file changed, 21 insertions(+) diff --git a/CHANGES.md b/CHANGES.md index dd315bc074..ed65cbc0b4 100644 --- a/CHANGES.md +++ b/CHANGES.md @@ -1,3 +1,24 @@ +# 24.11 + + * a new well-balanced method was added to the CTU PPM solver. This + does the characteristic projection only on the perturbed pressure + and then adds back in the hydrostatic pressure. It can be enabled + via `castro.ppm_well_balanced` (#2945) + + * fixed a bug in the div{U} calculation for artificial viscosity on + symmetry boundaries (#2983) + + * more development on 2D spherical geometry (#2973, #2975, #2981) + + * updates to the massive star plotting scripts (#2979) + + * add some new checks to prevent running unsupported combinations + of solvers (#2978) + + * documentation updates (#2977) + + * `flame_wave` can now be run in 1D (#2976) + # 24.10 * update initial model for `subchandra` when doing ASE NSE (#2970) From b38a018a1ee36964348a686a5dd8261a3056d050 Mon Sep 17 00:00:00 2001 From: Zhi Chen <62574124+zhichen3@users.noreply.github.com> Date: Fri, 8 Nov 2024 09:53:38 -0500 Subject: [PATCH 3/4] name variable to geom when its amrex::Geometry to avoid confusion (#2989) --- Source/driver/Derive.H | 70 ++++++++++++------------ Source/driver/Derive.cpp | 114 +++++++++++++++++++-------------------- 2 files changed, 92 insertions(+), 92 deletions(-) diff --git a/Source/driver/Derive.H b/Source/driver/Derive.H index 7dd302da16..0533c2d824 100644 --- a/Source/driver/Derive.H +++ b/Source/driver/Derive.H @@ -12,125 +12,125 @@ extern "C" void ca_derpres (const amrex::Box& bx, amrex::FArrayBox& derfab, int dcomp, int /*ncomp*/, - const amrex::FArrayBox& datfab, const amrex::Geometry& geomdata, + const amrex::FArrayBox& datfab, const amrex::Geometry& geom, amrex::Real /*time*/, const int* /*bcrec*/, int /*level*/); void ca_dereint1 (const amrex::Box& bx, amrex::FArrayBox& derfab, int dcomp, int /*ncomp*/, - const amrex::FArrayBox& datfab, const amrex::Geometry& geomdata, + const amrex::FArrayBox& datfab, const amrex::Geometry& geom, amrex::Real /*time*/, const int* /*bcrec*/, int /*level*/); void ca_dereint2 (const amrex::Box& bx, amrex::FArrayBox& derfab, int dcomp, int /*ncomp*/, - const amrex::FArrayBox& datfab, const amrex::Geometry& geomdata, + const amrex::FArrayBox& datfab, const amrex::Geometry& geom, amrex::Real /*time*/, const int* /*bcrec*/, int /*level*/); void ca_derlogden (const amrex::Box& bx, amrex::FArrayBox& derfab, int dcomp, int /*ncomp*/, - const amrex::FArrayBox& datfab, const amrex::Geometry& geomdata, + const amrex::FArrayBox& datfab, const amrex::Geometry& geom, amrex::Real /*time*/, const int* /*bcrec*/, int /*level*/); void ca_deruplusc (const amrex::Box& bx, amrex::FArrayBox& derfab, int dcomp, int /*ncomp*/, - const amrex::FArrayBox& datfab, const amrex::Geometry& geomdata, + const amrex::FArrayBox& datfab, const amrex::Geometry& geom, amrex::Real /*time*/, const int* /*bcrec*/, int /*level*/); void ca_deruminusc (const amrex::Box& bx, amrex::FArrayBox& derfab, int dcomp, int /*ncomp*/, - const amrex::FArrayBox& datfab, const amrex::Geometry& geomdata, + const amrex::FArrayBox& datfab, const amrex::Geometry& geom, amrex::Real /*time*/, const int* /*bcrec*/, int /*level*/); void ca_dersoundspeed (const amrex::Box& bx, amrex::FArrayBox& derfab, int dcomp, int /*ncomp*/, - const amrex::FArrayBox& datfab, const amrex::Geometry& geomdata, + const amrex::FArrayBox& datfab, const amrex::Geometry& geom, amrex::Real /*time*/, const int* /*bcrec*/, int /*level*/); void ca_dergamma1 (const amrex::Box& bx, amrex::FArrayBox& derfab, int dcomp, int /*ncomp*/, - const amrex::FArrayBox& datfab, const amrex::Geometry& geomdata, + const amrex::FArrayBox& datfab, const amrex::Geometry& geom, amrex::Real /*time*/, const int* /*bcrec*/, int /*level*/); void ca_dermachnumber (const amrex::Box& bx, amrex::FArrayBox& derfab, int dcomp, int /*ncomp*/, - const amrex::FArrayBox& datfab, const amrex::Geometry& geomdata, + const amrex::FArrayBox& datfab, const amrex::Geometry& geom, amrex::Real /*time*/, const int* /*bcrec*/, int /*level*/); void ca_derentropy (const amrex::Box& bx, amrex::FArrayBox& derfab, int dcomp, int /*ncomp*/, - const amrex::FArrayBox& datfab, const amrex::Geometry& geomdata, + const amrex::FArrayBox& datfab, const amrex::Geometry& geom, amrex::Real /*time*/, const int* /*bcrec*/, int /*level*/); #ifdef DIFFUSION void ca_dercond (const amrex::Box& bx, amrex::FArrayBox& derfab, int dcomp, int /*ncomp*/, - const amrex::FArrayBox& datfab, const amrex::Geometry& geomdata, + const amrex::FArrayBox& datfab, const amrex::Geometry& geom, amrex::Real /*time*/, const int* /*bcrec*/, int /*level*/); void ca_derdiffcoeff (const amrex::Box& bx, amrex::FArrayBox& derfab, int dcomp, int /*ncomp*/, - const amrex::FArrayBox& datfab, const amrex::Geometry& geomdata, + const amrex::FArrayBox& datfab, const amrex::Geometry& geom, amrex::Real /*time*/, const int* /*bcrec*/, int /*level*/); void ca_derdiffterm (const amrex::Box& bx, amrex::FArrayBox& derfab, int dcomp, int /*ncomp*/, - const amrex::FArrayBox& datfab, const amrex::Geometry& geomdata, + const amrex::FArrayBox& datfab, const amrex::Geometry& geom, amrex::Real /*time*/, const int* /*bcrec*/, int /*level*/); #endif void ca_derenuc (const amrex::Box& bx, amrex::FArrayBox& derfab, int dcomp, int /*ncomp*/, - const amrex::FArrayBox& datfab, const amrex::Geometry& geomdata, + const amrex::FArrayBox& datfab, const amrex::Geometry& geom, amrex::Real /*time*/, const int* /*bcrec*/, int /*level*/); void ca_derenuctimescale (const amrex::Box& bx, amrex::FArrayBox& derfab, int dcomp, int /*ncomp*/, - const amrex::FArrayBox& datfab, const amrex::Geometry& geomdata, + const amrex::FArrayBox& datfab, const amrex::Geometry& geom, amrex::Real /*time*/, const int* /*bcrec*/, int /*level*/); void ca_dervel (const amrex::Box& bx, amrex::FArrayBox& derfab, int dcomp, int /*ncomp*/, - const amrex::FArrayBox& datfab, const amrex::Geometry& geomdata, + const amrex::FArrayBox& datfab, const amrex::Geometry& geom, amrex::Real /*time*/, const int* /*bcrec*/, int /*level*/); void ca_dermagvel (const amrex::Box& bx, amrex::FArrayBox& derfab, int dcomp, int /*ncomp*/, - const amrex::FArrayBox& datfab, const amrex::Geometry& geomdata, + const amrex::FArrayBox& datfab, const amrex::Geometry& geom, amrex::Real /*time*/, const int* /*bcrec*/, int /*level*/); void ca_dermaggrav (const amrex::Box& bx, amrex::FArrayBox& derfab, int dcomp, int /*ncomp*/, - const amrex::FArrayBox& datfab, const amrex::Geometry& geomdata, + const amrex::FArrayBox& datfab, const amrex::Geometry& geom, amrex::Real /*time*/, const int* /*bcrec*/, int /*level*/); void ca_derradialvel (const amrex::Box& bx, amrex::FArrayBox& derfab, int dcomp, int /*ncomp*/, - const amrex::FArrayBox& datfab, const amrex::Geometry& geomdata, + const amrex::FArrayBox& datfab, const amrex::Geometry& geom, amrex::Real /*time*/, const int* /*bcrec*/, int /*level*/); void ca_dercircvel (const amrex::Box& bx, amrex::FArrayBox& derfab, int dcomp, int /*ncomp*/, - const amrex::FArrayBox& datfab, const amrex::Geometry& geomdata, + const amrex::FArrayBox& datfab, const amrex::Geometry& geom, amrex::Real /*time*/, const int* /*bcrec*/, int /*level*/); void ca_dermagmom (const amrex::Box& bx, amrex::FArrayBox& derfab, int dcomp, int /*ncomp*/, - const amrex::FArrayBox& datfab, const amrex::Geometry& geomdata, + const amrex::FArrayBox& datfab, const amrex::Geometry& geom, amrex::Real /*time*/, const int* /*bcrec*/, int /*level*/); void ca_derangmomx (const amrex::Box& bx, amrex::FArrayBox& derfab, int dcomp, int ncomp, - const amrex::FArrayBox& datfab, const amrex::Geometry& geomdata, + const amrex::FArrayBox& datfab, const amrex::Geometry& geom, amrex::Real time, const int* bcrec, int level); void ca_derangmomy (const amrex::Box& bx, amrex::FArrayBox& derfab, int dcomp, int ncomp, - const amrex::FArrayBox& datfab, const amrex::Geometry& geomdata, + const amrex::FArrayBox& datfab, const amrex::Geometry& geom, amrex::Real time, const int* bcrec, int level); void ca_derangmomz (const amrex::Box& bx, amrex::FArrayBox& derfab, int dcomp, int ncomp, - const amrex::FArrayBox& datfab, const amrex::Geometry& geomdata, + const amrex::FArrayBox& datfab, const amrex::Geometry& geom, amrex::Real time, const int* bcrec, int level); void ca_derkineng (const amrex::Box& bx, amrex::FArrayBox& derfab, int dcomp, int ncomp, - const amrex::FArrayBox& datfab, const amrex::Geometry& geomdata, + const amrex::FArrayBox& datfab, const amrex::Geometry& geom, amrex::Real time, const int* bcrec, int level); void ca_dernull @@ -144,53 +144,53 @@ extern "C" void ca_derspec (const amrex::Box& bx, amrex::FArrayBox& derfab, int dcomp, int /*ncomp*/, - const amrex::FArrayBox& datfab, const amrex::Geometry& geomdata, + const amrex::FArrayBox& datfab, const amrex::Geometry& geom, amrex::Real /*time*/, const int* /*bcrec*/, int /*level*/); void ca_derabar (const amrex::Box& bx, amrex::FArrayBox& derfab, int dcomp, int /*ncomp*/, - const amrex::FArrayBox& datfab, const amrex::Geometry& geomdata, + const amrex::FArrayBox& datfab, const amrex::Geometry& geom, amrex::Real /*time*/, const int* /*bcrec*/, int /*level*/); void ca_derye (const amrex::Box& bx, amrex::FArrayBox& derfab, int dcomp, int /*ncomp*/, - const amrex::FArrayBox& datfab, const amrex::Geometry& geomdata, + const amrex::FArrayBox& datfab, const amrex::Geometry& geom, amrex::Real /*time*/, const int* /*bcrec*/, int /*level*/); void ca_dermagvort (const amrex::Box& bx, amrex::FArrayBox& derfab, int dcomp, int /*ncomp*/, - const amrex::FArrayBox& datfab, const amrex::Geometry& geomdata, + const amrex::FArrayBox& datfab, const amrex::Geometry& geom, amrex::Real /*time*/, const int* /*bcrec*/, int /*level*/); void ca_derdivu (const amrex::Box& bx, amrex::FArrayBox& derfab, int dcomp, int /*ncomp*/, - const amrex::FArrayBox& datfab, const amrex::Geometry& geomdata, + const amrex::FArrayBox& datfab, const amrex::Geometry& geom, amrex::Real /*time*/, const int* /*bcrec*/, int /*level*/); void ca_derstate (const amrex::Box& bx, amrex::FArrayBox& derfab, int dcomp, int /*ncomp*/, - const amrex::FArrayBox& datfab, const amrex::Geometry& geomdata, + const amrex::FArrayBox& datfab, const amrex::Geometry& geom, amrex::Real /*time*/, const int* /*bcrec*/, int /*level*/); #ifdef MHD void ca_dermagcenx (const amrex::Box& bx, amrex::FArrayBox& derfab, int dcomp, int /*ncomp*/, - const amrex::FArrayBox& datfab, const amrex::Geometry& geomdata, + const amrex::FArrayBox& datfab, const amrex::Geometry& geom, amrex::Real /*time*/, const int* /*bcrec*/, int /*level*/); void ca_dermagceny (const amrex::Box& bx, amrex::FArrayBox& derfab, int dcomp, int /*ncomp*/, - const amrex::FArrayBox& datfab, const amrex::Geometry& geomdata, + const amrex::FArrayBox& datfab, const amrex::Geometry& geom, amrex::Real /*time*/, const int* /*bcrec*/, int /*level*/); void ca_dermagcenz (const amrex::Box& bx, amrex::FArrayBox& derfab, int dcomp, int /*ncomp*/, - const amrex::FArrayBox& datfab, const amrex::Geometry& geomdata, + const amrex::FArrayBox& datfab, const amrex::Geometry& geom, amrex::Real /*time*/, const int* /*bcrec*/, int /*level*/); void ca_derdivb (const amrex::Box& bx, amrex::FArrayBox& derfab, int dcomp, int /*ncomp*/, - const amrex::FArrayBox& datfab, const amrex::Geometry& geomdata, + const amrex::FArrayBox& datfab, const amrex::Geometry& geom, amrex::Real /*time*/, const int* /*bcrec*/, int /*level*/); #endif diff --git a/Source/driver/Derive.cpp b/Source/driver/Derive.cpp index 3064fc3dd2..63d31449d1 100644 --- a/Source/driver/Derive.cpp +++ b/Source/driver/Derive.cpp @@ -21,7 +21,7 @@ extern "C" // need to explicitly synchronize after GPU kernels. void ca_derpres(const Box& bx, FArrayBox& derfab, int /*dcomp*/, int /*ncomp*/, - const FArrayBox& datfab, const Geometry& /*geomdata*/, + const FArrayBox& datfab, const Geometry& /*geom*/, Real /*time*/, const int* /*bcrec*/, int /*level*/) { @@ -54,7 +54,7 @@ extern "C" } void ca_dereint1(const Box& bx, FArrayBox& derfab, int /*dcomp*/, int /*ncomp*/, - const FArrayBox& datfab, const Geometry& /*geomdata*/, + const FArrayBox& datfab, const Geometry& /*geom*/, Real /*time*/, const int* /*bcrec*/, int /*level*/) { @@ -76,7 +76,7 @@ extern "C" } void ca_dereint2(const Box& bx, FArrayBox& derfab, int /*dcomp*/, int /*ncomp*/, - const FArrayBox& datfab, const Geometry& /*geomdata*/, + const FArrayBox& datfab, const Geometry& /*geom*/, Real /*time*/, const int* /*bcrec*/, int /*level*/) { @@ -92,7 +92,7 @@ extern "C" } void ca_derlogden(const Box& bx, FArrayBox& derfab, int /*dcomp*/, int /*ncomp*/, - const FArrayBox& datfab, const Geometry& /*geomdata*/, + const FArrayBox& datfab, const Geometry& /*geom*/, Real /*time*/, const int* /*bcrec*/, int /*level*/) { @@ -107,7 +107,7 @@ extern "C" } void ca_deruplusc(const Box& bx, FArrayBox& derfab, int /*dcomp*/, int /*ncomp*/, - const FArrayBox& datfab, const Geometry& /*geomdata*/, + const FArrayBox& datfab, const Geometry& /*geom*/, Real /*time*/, const int* /*bcrec*/, int /*level*/) { @@ -142,7 +142,7 @@ extern "C" } void ca_deruminusc(const Box& bx, FArrayBox& derfab, int /*dcomp*/, int /*ncomp*/, - const FArrayBox& datfab, const Geometry& /*geomdata*/, + const FArrayBox& datfab, const Geometry& /*geom*/, Real /*time*/, const int* /*bcrec*/, int /*level*/) { @@ -177,7 +177,7 @@ extern "C" } void ca_dersoundspeed(const Box& bx, FArrayBox& derfab, int /*dcomp*/, int /*ncomp*/, - const FArrayBox& datfab, const Geometry& /*geomdata*/, + const FArrayBox& datfab, const Geometry& /*geom*/, Real /*time*/, const int* /*bcrec*/, int /*level*/) { @@ -213,7 +213,7 @@ extern "C" void ca_dergamma1(const Box& bx, FArrayBox& derfab, int /*dcomp*/, int /*ncomp*/, - const FArrayBox& datfab, const Geometry& /*geomdata*/, + const FArrayBox& datfab, const Geometry& /*geom*/, Real /*time*/, const int* /*bcrec*/, int /*level*/) { @@ -248,7 +248,7 @@ extern "C" } void ca_dermachnumber(const Box& bx, FArrayBox& derfab, int /*dcomp*/, int /*ncomp*/, - const FArrayBox& datfab, const Geometry& /*geomdata*/, + const FArrayBox& datfab, const Geometry& /*geom*/, Real /*time*/, const int* /*bcrec*/, int /*level*/) { @@ -286,7 +286,7 @@ extern "C" } void ca_derentropy(const Box& bx, FArrayBox& derfab, int /*dcomp*/, int /*ncomp*/, - const FArrayBox& datfab, const Geometry& /*geomdata*/, + const FArrayBox& datfab, const Geometry& /*geom*/, Real /*time*/, const int* /*bcrec*/, int /*level*/) { @@ -321,7 +321,7 @@ extern "C" #ifdef DIFFUSION void ca_dercond(const Box& bx, FArrayBox& derfab, int /*dcomp*/, int /*ncomp*/, - const FArrayBox& datfab, const Geometry& /*geomdata*/, + const FArrayBox& datfab, const Geometry& /*geom*/, Real /*time*/, const int* /*bcrec*/, int /*level*/) { @@ -333,7 +333,7 @@ extern "C" } void ca_derdiffcoeff(const Box& bx, FArrayBox& derfab, int /*dcomp*/, int /*ncomp*/, - const FArrayBox& datfab, const Geometry& /*geomdata*/, + const FArrayBox& datfab, const Geometry& /*geom*/, Real /*time*/, const int* /*bcrec*/, int /*level*/) { @@ -345,7 +345,7 @@ extern "C" } void ca_derdiffterm(const Box& bx, FArrayBox& derfab, int /*dcomp*/, int /*ncomp*/, - const FArrayBox& datfab, const Geometry& geomdata, + const FArrayBox& datfab, const Geometry& geom, Real /*time*/, const int* /*bcrec*/, int /*level*/) { @@ -360,9 +360,9 @@ extern "C" fill_temp_cond(obx, dat, coeff_arr); - auto dx = geomdata.CellSizeArray(); - auto problo = geomdata.ProbLoArray(); - const int coord_type = geomdata.Coord(); + auto dx = geom.CellSizeArray(); + auto problo = geom.ProbLoArray(); + const int coord_type = geom.Coord(); amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept @@ -432,14 +432,14 @@ extern "C" #ifdef REACTIONS void ca_derenuctimescale(const Box& bx, FArrayBox& derfab, int /*dcomp*/, int /*ncomp*/, - const FArrayBox& datfab, const Geometry& geomdata, + const FArrayBox& datfab, const Geometry& geom, Real /*time*/, const int* /*bcrec*/, int /*level*/) { auto const dat = datfab.array(); auto const der = derfab.array(); - auto dx = geomdata.CellSizeArray(); + auto dx = geom.CellSizeArray(); Real dd = 0.0_rt; #if AMREX_SPACEDIM == 1 @@ -494,7 +494,7 @@ extern "C" } void ca_derenuc(const Box& bx, FArrayBox& derfab, int /*dcomp*/, int /*ncomp*/, - const FArrayBox& datfab, const Geometry& /*geomdata*/, + const FArrayBox& datfab, const Geometry& /*geom*/, Real /*time*/, const int* /*bcrec*/, int /*level*/) { auto const dat = datfab.array(); @@ -512,7 +512,7 @@ extern "C" #endif void ca_dervel(const Box& bx, FArrayBox& derfab, int /*dcomp*/, int /*ncomp*/, - const FArrayBox& datfab, const Geometry& /*geomdata*/, + const FArrayBox& datfab, const Geometry& /*geom*/, Real /*time*/, const int* /*bcrec*/, int /*level*/) { @@ -528,7 +528,7 @@ extern "C" void ca_dermagvel(const Box& bx, FArrayBox& derfab, int /*dcomp*/, int /*ncomp*/, - const FArrayBox& datfab, const Geometry& /*geomdata*/, + const FArrayBox& datfab, const Geometry& /*geom*/, Real /*time*/, const int* /*bcrec*/, int /*level*/) { @@ -549,7 +549,7 @@ extern "C" void ca_dermaggrav(const Box& bx, FArrayBox& derfab, int /*dcomp*/, int /*ncomp*/, - const FArrayBox& datfab, const Geometry& /*geomdata*/, + const FArrayBox& datfab, const Geometry& /*geom*/, Real /*time*/, const int* /*bcrec*/, int /*level*/) { @@ -568,7 +568,7 @@ extern "C" } void ca_derradialvel(const Box& bx, FArrayBox& derfab, int /*dcomp*/, int /*ncomp*/, - const FArrayBox& datfab, const Geometry& geomdata, + const FArrayBox& datfab, const Geometry& geom, Real /*time*/, const int* /*bcrec*/, int /*level*/) { @@ -577,8 +577,8 @@ extern "C" auto const dat = datfab.array(); auto const der = derfab.array(); - auto dx = geomdata.CellSizeArray(); - auto problo = geomdata.ProbLoArray(); + auto dx = geom.CellSizeArray(); + auto problo = geom.ProbLoArray(); amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept @@ -623,7 +623,7 @@ extern "C" void ca_dercircvel(const Box& bx, FArrayBox& derfab, int /*dcomp*/, int /*ncomp*/, - const FArrayBox& datfab, const Geometry& geomdata, + const FArrayBox& datfab, const Geometry& geom, Real /*time*/, const int* /*bcrec*/, int /*level*/) { @@ -632,8 +632,8 @@ extern "C" auto const dat = datfab.array(); auto const der = derfab.array(); - auto dx = geomdata.CellSizeArray(); - auto problo = geomdata.ProbLoArray(); + auto dx = geom.CellSizeArray(); + auto problo = geom.ProbLoArray(); amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept @@ -688,7 +688,7 @@ extern "C" void ca_dermagmom(const Box& bx, FArrayBox& derfab, int /*dcomp*/, int /*ncomp*/, - const FArrayBox& datfab, const Geometry& /*geomdata*/, + const FArrayBox& datfab, const Geometry& /*geom*/, Real /*time*/, const int* /*bcrec*/, int /*level*/) { @@ -707,13 +707,13 @@ extern "C" } void ca_derangmomx (const Box& bx, FArrayBox& derfab, int /*dcomp*/, int /*ncomp*/, - const FArrayBox& datfab, const Geometry& geomdata, + const FArrayBox& datfab, const Geometry& geom, Real /*time*/, const int* /*bcrec*/, int /*level*/) { int idir = 0; - auto dx = geomdata.CellSizeArray(); - auto problo = geomdata.ProbLoArray(); + auto dx = geom.CellSizeArray(); + auto problo = geom.ProbLoArray(); auto const dat = datfab.array(); auto const L = derfab.array(); @@ -724,7 +724,7 @@ extern "C" Real loc[3]; //loc calculated like sum_utils.cpp - //This might be equivalent and more modular: position(i, j, k, geomdata, loc); + //This might be equivalent and more modular: position(i, j, k, geom, loc); loc[0] = problo[0] + (0.5_rt + i) * dx[0]; #if AMREX_SPACEDIM >= 2 @@ -759,13 +759,13 @@ extern "C" } void ca_derangmomy (const Box& bx, FArrayBox& derfab, int /*dcomp*/, int /*ncomp*/, - const FArrayBox& datfab, const Geometry& geomdata, + const FArrayBox& datfab, const Geometry& geom, Real /*time*/, const int* /*bcrec*/, int /*level*/) { int idir = 1; - auto dx = geomdata.CellSizeArray(); - auto problo = geomdata.ProbLoArray(); + auto dx = geom.CellSizeArray(); + auto problo = geom.ProbLoArray(); auto const dat = datfab.array(); auto const L = derfab.array(); @@ -807,13 +807,13 @@ extern "C" } void ca_derangmomz (const Box& bx, FArrayBox& derfab, int /*dcomp*/, int /*ncomp*/, - const FArrayBox& datfab, const Geometry& geomdata, + const FArrayBox& datfab, const Geometry& geom, Real /*time*/, const int* /*bcrec*/, int /*level*/) { int idir = 2; - auto dx = geomdata.CellSizeArray(); - auto problo = geomdata.ProbLoArray(); + auto dx = geom.CellSizeArray(); + auto problo = geom.ProbLoArray(); auto const dat = datfab.array(); auto const L = derfab.array(); @@ -856,7 +856,7 @@ extern "C" } void ca_derkineng (const Box& bx, FArrayBox& derfab, int /*dcomp*/, int /*ncomp*/, - const FArrayBox& datfab, const Geometry& /*geomdata*/, + const FArrayBox& datfab, const Geometry& /*geom*/, Real /*time*/, const int* /*bcrec*/, int /*level*/) { @@ -887,7 +887,7 @@ extern "C" } void ca_derspec(const Box& bx, FArrayBox& derfab, int /*dcomp*/, int /*ncomp*/, - const FArrayBox& datfab, const Geometry& /*geomdata*/, + const FArrayBox& datfab, const Geometry& /*geom*/, Real /*time*/, const int* /*bcrec*/, int /*level*/) { @@ -903,7 +903,7 @@ extern "C" void ca_derye(const Box& bx, FArrayBox& derfab, int /*dcomp*/, int /*ncomp*/, - const FArrayBox& datfab, const Geometry& /*geomdata*/, + const FArrayBox& datfab, const Geometry& /*geom*/, Real /*time*/, const int* /*bcrec*/, int /*level*/) { @@ -925,7 +925,7 @@ extern "C" } void ca_derabar(const Box& bx, FArrayBox& derfab, int /*dcomp*/, int /*ncomp*/, - const FArrayBox& datfab, const Geometry& /*geomdata*/, + const FArrayBox& datfab, const Geometry& /*geom*/, Real /*time*/, const int* /*bcrec*/, int /*level*/) { @@ -947,19 +947,19 @@ extern "C" } void ca_dermagvort(const Box& bx, FArrayBox& derfab, int /*dcomp*/, int /*ncomp*/, - const FArrayBox& datfab, const Geometry& geomdata, + const FArrayBox& datfab, const Geometry& geom, Real /*time*/, const int* /*bcrec*/, int /*level*/) { auto const dat = datfab.array(); auto const der = derfab.array(); - auto dx = geomdata.CellSizeArray(); + auto dx = geom.CellSizeArray(); - const int coord_type = geomdata.Coord(); + const int coord_type = geom.Coord(); #if AMREX_SPACEDIM == 2 - auto problo = geomdata.ProbLoArray(); + auto problo = geom.ProbLoArray(); #endif amrex::ParallelFor(bx, @@ -1046,18 +1046,18 @@ extern "C" } void ca_derdivu(const Box& bx, FArrayBox& derfab, int /*dcomp*/, int /*ncomp*/, - const FArrayBox& datfab, const Geometry& geomdata, + const FArrayBox& datfab, const Geometry& geom, Real /*time*/, const int* /*bcrec*/, int /*level*/) { auto const dat = datfab.array(); auto const der = derfab.array(); - auto dx = geomdata.CellSizeArray(); + auto dx = geom.CellSizeArray(); - auto problo = geomdata.ProbLoArray(); + auto problo = geom.ProbLoArray(); - const int coord_type = geomdata.Coord(); + const int coord_type = geom.Coord(); amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept @@ -1112,7 +1112,7 @@ extern "C" } void ca_derstate(const Box& bx, FArrayBox& derfab, int /*dcomp*/, int /*ncomp*/, - const FArrayBox& datfab, const Geometry& /*geomdata*/, + const FArrayBox& datfab, const Geometry& /*geom*/, Real /*time*/, const int* /*bcrec*/, int /*level*/) { @@ -1138,7 +1138,7 @@ extern "C" #ifdef MHD void ca_dermagcenx(const Box& bx, FArrayBox& derfab, int /*dcomp*/, int /*ncomp*/, - const FArrayBox& datfab, const Geometry& /*geomdata*/, + const FArrayBox& datfab, const Geometry& /*geom*/, Real /*time*/, const int* /*bcrec*/, int /*level*/) { @@ -1154,7 +1154,7 @@ extern "C" } void ca_dermagceny(const Box& bx, FArrayBox& derfab, int /*dcomp*/, int /*ncomp*/, - const FArrayBox& datfab, const Geometry& /*geomdata*/, + const FArrayBox& datfab, const Geometry& /*geom*/, Real /*time*/, const int* /*bcrec*/, int /*level*/) { @@ -1170,7 +1170,7 @@ extern "C" } void ca_dermagcenz(const Box& bx, FArrayBox& derfab, int /*dcomp*/, int /*ncomp*/, - const FArrayBox& datfab, const Geometry& /*geomdata*/, + const FArrayBox& datfab, const Geometry& /*geom*/, Real /*time*/, const int* /*bcrec*/, int /*level*/) { @@ -1186,14 +1186,14 @@ extern "C" } void ca_derdivb(const Box& bx, FArrayBox& derfab, int /*dcomp*/, int /*ncomp*/, - const FArrayBox& datfab, const Geometry& geomdata, + const FArrayBox& datfab, const Geometry& geom, Real /*time*/, const int* /*bcrec*/, int /*level*/) { auto const dat = datfab.array(); auto const der = derfab.array(); - auto dx = geomdata.CellSizeArray(); + auto dx = geom.CellSizeArray(); amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept From 161c875ae969b327d0aef545145278f8a61ba1e0 Mon Sep 17 00:00:00 2001 From: Zhi Chen <62574124+zhichen3@users.noreply.github.com> Date: Wed, 13 Nov 2024 11:56:18 -0500 Subject: [PATCH 4/4] add a distance function (#2990) Adds a separate function to calculate distance, mainly because the old way gets the 2d spherical case wrong. --- Source/driver/Castro.cpp | 17 ++++++----- Source/driver/Castro_util.H | 22 ++++++++++++++ Source/driver/Derive.cpp | 49 ++++++++++++++----------------- Source/gravity/Gravity.cpp | 13 ++++---- Source/reactions/Castro_react.cpp | 10 ++++--- Source/sources/Castro_sponge.cpp | 3 +- 6 files changed, 69 insertions(+), 45 deletions(-) diff --git a/Source/driver/Castro.cpp b/Source/driver/Castro.cpp index 6149e7f988..2a47ac55dd 100644 --- a/Source/driver/Castro.cpp +++ b/Source/driver/Castro.cpp @@ -3641,19 +3641,17 @@ Castro::apply_tagging_restrictions(TagBoxArray& tags, [[maybe_unused]] Real time const Real* probhi = geomdata.ProbHi(); const Real* dx = geomdata.CellSize(); - Real loc[3] = {0.0}; + GpuArray loc = {0.0}; - loc[0] = problo[0] + (static_cast(i) + 0.5_rt) * dx[0]; + loc[0] = problo[0] + (static_cast(i) + 0.5_rt) * dx[0] - problem::center[0]; #if AMREX_SPACEDIM >= 2 - loc[1] = problo[1] + (static_cast(j) + 0.5_rt) * dx[1]; + loc[1] = problo[1] + (static_cast(j) + 0.5_rt) * dx[1] - problem::center[1]; #endif #if AMREX_SPACEDIM == 3 - loc[2] = problo[2] + (static_cast(k) + 0.5_rt) * dx[2]; + loc[2] = problo[2] + (static_cast(k) + 0.5_rt) * dx[2] - problem::center[2]; #endif - Real r = std::sqrt((loc[0] - problem::center[0]) * (loc[0] - problem::center[0]) + - (loc[1] - problem::center[1]) * (loc[1] - problem::center[1]) + - (loc[2] - problem::center[2]) * (loc[2] - problem::center[2])); + Real r = distance(geomdata, loc); Real max_dist_lo = 0.0; Real max_dist_hi = 0.0; @@ -4357,9 +4355,12 @@ Castro::define_new_center(const MultiFab& S, Real time) // Now broadcast to everyone else. ParallelDescriptor::Bcast(&problem::center[0], AMREX_SPACEDIM, owner); - // Make sure if R-Z that center stays exactly on axis + // Make sure if R-Z and SPHERICAL that center stays exactly on axis if ( Geom().IsRZ() ) { problem::center[0] = 0; + } else if ( Geom().IsSPHERICAL() ) { + problem::center[0] = 0; + problem::center[1] = 0; } } diff --git a/Source/driver/Castro_util.H b/Source/driver/Castro_util.H index 016215aa32..c6162b0190 100644 --- a/Source/driver/Castro_util.H +++ b/Source/driver/Castro_util.H @@ -156,6 +156,28 @@ void position(int i, int j, int k, } + +AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE +Real distance(GeometryData const& geomdata, GpuArray& loc) +{ + // Returns the distance from the center provided with loc, the position. + // loc is the form of [x, y, z,] in Cartesian, [r, z, phi] in cylindrical + // and [r, theta, phi] in spherical + + const auto coord = geomdata.Coord(); + + if (coord == CoordSys::cartesian) { + return std::sqrt(loc[0]*loc[0] + loc[1]*loc[1] + loc[2]*loc[2]); + } + + if (coord == CoordSys::RZ) { + return std::sqrt(loc[0]*loc[0] + loc[1]*loc[1]); + } + + return std::abs(loc[0]); +} + + namespace geometry_util { diff --git a/Source/driver/Derive.cpp b/Source/driver/Derive.cpp index 63d31449d1..78a1dfcc35 100644 --- a/Source/driver/Derive.cpp +++ b/Source/driver/Derive.cpp @@ -579,21 +579,18 @@ extern "C" auto dx = geom.CellSizeArray(); auto problo = geom.ProbLoArray(); + auto geomdata = geom.data(); amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept { - - Real x = problo[0] + (static_cast(i) + 0.5_rt) * dx[0] - problem::center[0]; + GpuArray loc = {0.0}; + loc[0] = problo[0] + (static_cast(i) + 0.5_rt) * dx[0] - problem::center[0]; #if AMREX_SPACEDIM >= 2 - Real y = problo[1] + (static_cast(j) + 0.5_rt) * dx[1] - problem::center[1]; -#else - Real y = 0.0_rt; + loc[1] = problo[1] + (static_cast(j) + 0.5_rt) * dx[1] - problem::center[1]; #endif #if AMREX_SPACEDIM == 3 - Real z = problo[2] + (static_cast(k) + 0.5_rt) * dx[2] - problem::center[2]; -#else - Real z = 0.0_rt; + loc[2] = problo[2] + (static_cast(k) + 0.5_rt) * dx[2] - problem::center[2]; #endif if (domain_is_plane_parallel) { @@ -607,15 +604,15 @@ extern "C" // where e_r and e_phi are the cylindrical unit vectors // we need the distance in the x-y plane from the origin - Real r = std::sqrt(x*x + y*y); - der(i,j,k,0) = (dat(i,j,k,1)*x + dat(i,j,k,2)*y) / (dat(i,j,k,0)*r); + Real r = std::sqrt(loc[0]*loc[0] + loc[1]*loc[1]); + der(i,j,k,0) = (dat(i,j,k,1)*loc[0] + dat(i,j,k,2)*loc[1]) / (dat(i,j,k,0)*r); #endif } else { - Real r = std::sqrt(x*x + y*y + z*z); + Real r = distance(geomdata, loc); - der(i,j,k,0) = (dat(i,j,k,1)*x + - dat(i,j,k,2)*y + - dat(i,j,k,3)*z) / ( dat(i,j,k,0)*r ); + der(i,j,k,0) = (dat(i,j,k,1)*loc[0] + + dat(i,j,k,2)*loc[1] + + dat(i,j,k,3)*loc[2]) / ( dat(i,j,k,0)*r ); } }); @@ -634,21 +631,19 @@ extern "C" auto dx = geom.CellSizeArray(); auto problo = geom.ProbLoArray(); + auto geomdata = geom.data(); amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept { - Real x = problo[0] + (static_cast(i) + 0.5_rt) * dx[0] - problem::center[0]; + GpuArray loc = {0.0}; + loc[0] = problo[0] + (static_cast(i) + 0.5_rt) * dx[0] - problem::center[0]; #if AMREX_SPACEDIM >= 2 - Real y = problo[1] + (static_cast(j) + 0.5_rt) * dx[1] - problem::center[1]; -#else - Real y = 0.0_rt; + loc[1] = problo[1] + (static_cast(j) + 0.5_rt) * dx[1] - problem::center[1]; #endif #if AMREX_SPACEDIM == 3 - Real z = problo[2] + (static_cast(k) + 0.5_rt) * dx[2] - problem::center[2]; -#else - Real z = 0.0_rt; + loc[2] = problo[2] + (static_cast(k) + 0.5_rt) * dx[2] - problem::center[2]; #endif if (domain_is_plane_parallel) { @@ -662,11 +657,11 @@ extern "C" // where e_r and e_phi are the cylindrical unit vectors // we need the distance in the x-y plane from the origin - Real r = std::sqrt(x*x + y*y); - der(i,j,k,0) = (-dat(i,j,k,1)*y + dat(i,j,k,2)*x) / (dat(i,j,k,0)*r); + Real r = std::sqrt(loc[0]*loc[0] + loc[1]*loc[1]); + der(i,j,k,0) = (-dat(i,j,k,1)*loc[1] + dat(i,j,k,2)*loc[0]) / (dat(i,j,k,0)*r); #endif } else { - Real r = std::sqrt(x*x + y*y + z*z); + Real r = distance(geomdata, loc); // we really mean just the velocity component that is // perpendicular to radial, and in general 3-d (e.g. a @@ -676,9 +671,9 @@ extern "C" dat(i,j,k,2)*dat(i,j,k,2) + dat(i,j,k,3)*dat(i,j,k,3))/(dat(i,j,k,0)*dat(i,j,k,0)); - Real vr = (dat(i,j,k,1)*x + - dat(i,j,k,2)*y + - dat(i,j,k,3)*z) / ( dat(i,j,k,0)*r ); + Real vr = (dat(i,j,k,1)*loc[0] + + dat(i,j,k,2)*loc[1] + + dat(i,j,k,3)*loc[2]) / ( dat(i,j,k,0)*r ); der(i,j,k,0) = std::sqrt(amrex::max(vtot2 - vr*vr, 0.0_rt)); } diff --git a/Source/gravity/Gravity.cpp b/Source/gravity/Gravity.cpp index 82003be86f..cfd40f65f0 100644 --- a/Source/gravity/Gravity.cpp +++ b/Source/gravity/Gravity.cpp @@ -1334,6 +1334,7 @@ Gravity::interpolate_monopole_grav(int level, RealVector& radial_grav, MultiFab& const Real dr = dx[0] / static_cast(gravity::drdxfac); const auto problo = geom.ProbLoArray(); + const auto geomdata = geom.data(); #ifdef _OPENMP #pragma omp parallel @@ -1367,7 +1368,7 @@ Gravity::interpolate_monopole_grav(int level, RealVector& radial_grav, MultiFab& loc[2] = 0.0_rt; #endif - Real r = std::sqrt(loc[0] * loc[0] + loc[1] * loc[1] + loc[2] * loc[2]); + Real r = distance(geomdata, loc); int index = static_cast(r / dr); @@ -1457,6 +1458,7 @@ Gravity::compute_radial_mass(const Box& bx, Real drinv = 1.0_rt / dr; const int coord_type = geom.Coord(); + const auto geomdata = geom.data(); AMREX_ALWAYS_ASSERT(coord_type >= 0 && coord_type <= 2); @@ -1497,16 +1499,17 @@ Gravity::compute_radial_mass(const Box& bx, amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept { - Real xc = problo[0] + (static_cast(i) + 0.5_rt) * dx[0] - problem::center[0]; + GpuArray loc; + loc[0] = problo[0] + (static_cast(i) + 0.5_rt) * dx[0] - problem::center[0]; Real lo_i = problo[0] + static_cast(i) * dx[0] - problem::center[0]; - Real yc = problo[1] + (static_cast(j) + 0.5_rt) * dx[1] - problem::center[1]; + loc[1]= problo[1] + (static_cast(j) + 0.5_rt) * dx[1] - problem::center[1]; Real lo_j = problo[1] + static_cast(j) * dx[1] - problem::center[1]; - Real zc = problo[2] + (static_cast(k) + 0.5_rt) * dx[2] - problem::center[2]; + loc[2]= problo[2] + (static_cast(k) + 0.5_rt) * dx[2] - problem::center[2]; Real lo_k = problo[2] + static_cast(k) * dx[2] - problem::center[2]; - Real r = std::sqrt(xc * xc + yc * yc + zc * zc); + Real r = distance(geomdata, loc); int index = static_cast(r * drinv); // We may be coming in here with a masked out zone (in a zone on a coarse diff --git a/Source/reactions/Castro_react.cpp b/Source/reactions/Castro_react.cpp index b7a422f1dd..9927bc4f2a 100644 --- a/Source/reactions/Castro_react.cpp +++ b/Source/reactions/Castro_react.cpp @@ -209,6 +209,7 @@ Castro::react_state(MultiFab& s, MultiFab& r, Real time, Real dt, const int stra const auto dx = geom.CellSizeArray(); #ifdef MODEL_PARSER const auto problo = geom.ProbLoArray(); + const auto geomdata = geom.data(); #endif #if defined(AMREX_USE_GPU) @@ -270,7 +271,7 @@ Castro::react_state(MultiFab& s, MultiFab& r, Real time, Real dt, const int stra #ifdef MODEL_PARSER if (drive_initial_convection) { - Real rr[3] = {0.0_rt}; + GpuArray rr = {0.0_rt}; rr[0] = problo[0] + dx[0] * (static_cast(i) + 0.5_rt) - problem::center[0]; #if AMREX_SPACEDIM >= 2 @@ -285,7 +286,7 @@ Castro::react_state(MultiFab& s, MultiFab& r, Real time, Real dt, const int stra if (domain_is_plane_parallel) { dist = rr[AMREX_SPACEDIM-1]; } else { - dist = std::sqrt(rr[0] * rr[0] + rr[1] * rr[1] + rr[2] * rr[2]); + dist = distance(geomdata, rr); } burn_state.T_fixed = interpolate(dist, model::itemp); @@ -565,6 +566,7 @@ Castro::react_state(Real time, Real dt) const auto dx = geom.CellSizeArray(); const auto problo = geom.ProbLoArray(); + const auto geomdata = geom.data(); #if defined(AMREX_USE_GPU) ParallelFor(bx, [=] AMREX_GPU_DEVICE (int i, int j, int k) @@ -618,7 +620,7 @@ Castro::react_state(Real time, Real dt) #ifdef MODEL_PARSER if (drive_initial_convection) { - Real rr[3] = {0.0_rt}; + GpuArray rr = {0.0_rt}; rr[0] = problo[0] + dx[0] * (static_cast(i) + 0.5_rt) - problem::center[0]; #if AMREX_SPACEDIM >= 2 @@ -633,7 +635,7 @@ Castro::react_state(Real time, Real dt) if (domain_is_plane_parallel) { dist = rr[AMREX_SPACEDIM-1]; } else { - dist = std::sqrt(rr[0] * rr[0] + rr[1] * rr[1] + rr[2] * rr[2]); + dist = distance(geomdata, rr); } burn_state.T_fixed = interpolate(dist, model::itemp); diff --git a/Source/sources/Castro_sponge.cpp b/Source/sources/Castro_sponge.cpp index 2109d63374..06781f4403 100644 --- a/Source/sources/Castro_sponge.cpp +++ b/Source/sources/Castro_sponge.cpp @@ -81,6 +81,7 @@ Castro::apply_sponge(const Box& bx, auto dx = geom.CellSizeArray(); auto problo = geom.ProbLoArray(); + auto geomdata = geom.data(); amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept @@ -124,7 +125,7 @@ Castro::apply_sponge(const Box& bx, Real sponge_factor = 0.0_rt; if (sponge_lower_radius >= 0.0_rt && sponge_upper_radius > sponge_lower_radius) { - Real rad = std::sqrt(r[0]*r[0] + r[1]*r[1] + r[2]*r[2]); + Real rad = distance(geomdata, r); if (rad < sponge_lower_radius) { sponge_factor = sponge_lower_factor;