From 4eabeab9f2f1183515758abc205dbc411c748bb8 Mon Sep 17 00:00:00 2001 From: Michael Zingale Date: Tue, 30 May 2023 14:22:22 -0400 Subject: [PATCH 01/11] add geometric source terms to the PPM tracing --- Source/hydro/reconstruction.H | 17 +++++++++++++++++ Source/hydro/trace_ppm.cpp | 28 +++++++++++++++++++++++++--- 2 files changed, 42 insertions(+), 3 deletions(-) diff --git a/Source/hydro/reconstruction.H b/Source/hydro/reconstruction.H index 6ee405a39f..1cd6a03c38 100644 --- a/Source/hydro/reconstruction.H +++ b/Source/hydro/reconstruction.H @@ -78,5 +78,22 @@ load_passive_stencil(Array4 const& U_arr, Array4 const& } +AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE +void +add_geometric_rho_source(Array4 const& q_arr, + Array4 const& dloga, + const int i, const int j, const int k, + Real* s) { + + // note: this is assumed to be working only in the x-direction + + s[im2] += -dloga(i-2,j,k) * q_arr(i-2,j,k,QRHO) * q_arr(i-2,j,k,QU); + s[im1] += -dloga(i-1,j,k) * q_arr(i-1,j,k,QRHO) * q_arr(i-1,j,k,QU); + s[i0] += -dloga(i,j,k) * q_arr(i,j,k,QRHO) * q_arr(i,j,k,QU); + s[ip1] += -dloga(i+1,j,k) * q_arr(i+1,j,k,QRHO) * q_arr(i+1,j,k,QU); + s[ip2] += -dloga(i+2,j,k) * q_arr(i+2,j,k,QRHO) * q_arr(i+2,j,k,QU); +} + + #endif diff --git a/Source/hydro/trace_ppm.cpp b/Source/hydro/trace_ppm.cpp index 02c819c9fa..f3017c0de9 100644 --- a/Source/hydro/trace_ppm.cpp +++ b/Source/hydro/trace_ppm.cpp @@ -56,6 +56,7 @@ Castro::trace_ppm(const Box& bx, const auto dx = geom.CellSizeArray(); + const int coord = geom.Coord(); Real hdt = 0.5_rt * dt; Real dtdx = dt / dx[idir]; @@ -81,6 +82,12 @@ Castro::trace_ppm(const Box& bx, for (int n = 0; n < NQSRC; n++) { do_source_trace[n] = 0; + // geometric source terms in r-direction need tracing + if (coord > 0 && idir == 0 && (n == QRHO || n == QPRES || n == QRHOE)) { + do_source_trace[n] = 1; + continue; + } + for (int k = lo[2]-2*dg2; k <= hi[2]+2*dg2; k++) { for (int j = lo[1]-2*dg1; j <= hi[1]+2*dg1; j++) { for (int i = lo[0]-2; i <= hi[0]+2; i++) { @@ -278,11 +285,18 @@ Castro::trace_ppm(const Box& bx, #ifndef AMREX_USE_GPU do_trace = do_source_trace[QRHO]; #else - do_trace = check_trace_source(srcQ, idir, i, j, k, QRHO); + if (idir == 0 && coord > 0) { + do_trace = 1; + } else { + do_trace = check_trace_source(srcQ, idir, i, j, k, QRHO); + } #endif if (do_trace) { load_stencil(srcQ, idir, i, j, k, QRHO, s); + if (idir == 0 && coord > 0) { + add_geometric_rho_source(q_arr, dloga, i, j, k, s); + } ppm_reconstruct(s, flat, sm, sp); ppm_int_profile(sm, sp, s[i0], un, cc, dtdx, Ip_src_rho, Im_src_rho); } @@ -315,7 +329,11 @@ Castro::trace_ppm(const Box& bx, #ifndef AMREX_USE_GPU do_trace = do_source_trace[QPRES]; #else - do_trace = check_trace_source(srcQ, idir, i, j, k, QPRES); + if (idir == 0 && coord > 0) { + do_trace = 1; + } else { + do_trace = check_trace_source(srcQ, idir, i, j, k, QPRES); + } #endif if (do_trace) { @@ -332,7 +350,11 @@ Castro::trace_ppm(const Box& bx, #ifndef AMREX_USE_GPU do_trace = do_source_trace[QREINT]; #else - do_trace = check_trace_source(srcQ, idir, i, j, k, QREINT); + if (idir == 0 && coord > 0) { + do_trace = 1; + } else { + do_trace = check_trace_source(srcQ, idir, i, j, k, QREINT); + } #endif if (do_trace) { From 2a44c103b6006a0223de12c1cbf31629bc439801 Mon Sep 17 00:00:00 2001 From: Michael Zingale Date: Tue, 30 May 2023 15:05:46 -0400 Subject: [PATCH 02/11] finish implementation --- Source/hydro/reconstruction.H | 39 +++++++++++++++++++++++++++++++++++ Source/hydro/trace_ppm.cpp | 37 +++++++-------------------------- 2 files changed, 46 insertions(+), 30 deletions(-) diff --git a/Source/hydro/reconstruction.H b/Source/hydro/reconstruction.H index 1cd6a03c38..65d6fc31dc 100644 --- a/Source/hydro/reconstruction.H +++ b/Source/hydro/reconstruction.H @@ -85,6 +85,8 @@ add_geometric_rho_source(Array4 const& q_arr, const int i, const int j, const int k, Real* s) { + using namespace reconstruction; + // note: this is assumed to be working only in the x-direction s[im2] += -dloga(i-2,j,k) * q_arr(i-2,j,k,QRHO) * q_arr(i-2,j,k,QU); @@ -94,6 +96,43 @@ add_geometric_rho_source(Array4 const& q_arr, s[ip2] += -dloga(i+2,j,k) * q_arr(i+2,j,k,QRHO) * q_arr(i+2,j,k,QU); } +AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE +void +add_geometric_rhoe_source(Array4 const& q_arr, + Array4 const& dloga, + const int i, const int j, const int k, + Real* s) { + + using namespace reconstruction; + + // note: this is assumed to be working only in the x-direction + + s[im2] += -dloga(i-2,j,k) * (q_arr(i-2,j,k,QREINT) + q_arr(i-2,j,k,QPRES)) * q_arr(i-2,j,k,QU); + s[im1] += -dloga(i-1,j,k) * (q_arr(i-1,j,k,QREINT) + q_arr(i-1,j,k,QPRES)) * q_arr(i-1,j,k,QU); + s[i0] += -dloga(i,j,k) * (q_arr(i,j,k,QREINT) + q_arr(i,j,k,QPRES)) * q_arr(i,j,k,QU); + s[ip1] += -dloga(i+1,j,k) * (q_arr(i+1,j,k,QREINT) + q_arr(i+1,j,k,QPRES)) * q_arr(i+1,j,k,QU); + s[ip2] += -dloga(i+2,j,k) * (q_arr(i+2,j,k,QREINT) + q_arr(i+2,j,k,QPRES)) * q_arr(i+2,j,k,QU); +} + +AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE +void +add_geometric_p_source(Array4 const& q_arr, + Array4 const& qaux_arr, + Array4 const& dloga, + const int i, const int j, const int k, + Real* s) { + + using namespace reconstruction; + + // note: this is assumed to be working only in the x-direction + + s[im2] += -dloga(i-2,j,k) * q_arr(i-2,j,k,QRHO) * qaux_arr(i-2,j,k,QC) * q_arr(i-2,j,k,QU); + s[im1] += -dloga(i-1,j,k) * q_arr(i-1,j,k,QRHO) * qaux_arr(i-1,j,k,QC) * q_arr(i-1,j,k,QU); + s[i0] += -dloga(i,j,k) * q_arr(i,j,k,QRHO) * qaux_arr(i,j,k,QC) * q_arr(i,j,k,QU); + s[ip1] += -dloga(i+1,j,k) * q_arr(i+1,j,k,QRHO) * qaux_arr(i+1,j,k,QC) * q_arr(i+1,j,k,QU); + s[ip2] += -dloga(i+2,j,k) * q_arr(i+2,j,k,QRHO) * qaux_arr(i+2,j,k,QC) * q_arr(i+2,j,k,QU); +} + #endif diff --git a/Source/hydro/trace_ppm.cpp b/Source/hydro/trace_ppm.cpp index f3017c0de9..3f157305a9 100644 --- a/Source/hydro/trace_ppm.cpp +++ b/Source/hydro/trace_ppm.cpp @@ -83,7 +83,7 @@ Castro::trace_ppm(const Box& bx, do_source_trace[n] = 0; // geometric source terms in r-direction need tracing - if (coord > 0 && idir == 0 && (n == QRHO || n == QPRES || n == QRHOE)) { + if (coord > 0 && idir == 0 && (n == QRHO || n == QPRES || n == QREINT)) { do_source_trace[n] = 1; continue; } @@ -338,6 +338,9 @@ Castro::trace_ppm(const Box& bx, if (do_trace) { load_stencil(srcQ, idir, i, j, k, QPRES, s); + if (idir == 0 && coord > 0) { + add_geometric_p_source(q_arr, qaux_arr, dloga, i, j, k, s); + } ppm_reconstruct(s, flat, sm, sp); ppm_int_profile(sm, sp, s[i0], un, cc, dtdx, Ip_src_p, Im_src_p); } @@ -359,6 +362,9 @@ Castro::trace_ppm(const Box& bx, if (do_trace) { load_stencil(srcQ, idir, i, j, k, QREINT, s); + if (idir == 0 && coord > 0) { + add_geometric_rhoe_source(q_arr, dloga, i, j, k, s); + } ppm_reconstruct(s, flat, sm, sp); ppm_int_profile(sm, sp, s[i0], un, cc, dtdx, Ip_src_rhoe, Im_src_rhoe); } @@ -625,35 +631,6 @@ Castro::trace_ppm(const Box& bx, } - // geometry source terms -#if (AMREX_SPACEDIM < 3) - // these only apply for x states (idir = 0) - if (idir == 0 && dloga(i,j,k) != 0.0_rt) { - Real rho = q_arr(i,j,k,QRHO); - - Real courn = dt/dx[0]*(cc+std::abs(un)); - Real eta = (1.0_rt - courn)/(cc*dt*std::abs(dloga(i,j,k))); - Real dlogatmp = amrex::min(eta, 1.0_rt)*dloga(i,j,k); - Real sourcr = -0.5_rt*dt*rho*dlogatmp*un; - Real sourcp = sourcr*csq; - Real source = sourcp*((q_arr(i,j,k,QPRES) + q_arr(i,j,k,QREINT))/rho)/csq; - - if (i <= vhi[0]) { - qm(i+1,j,k,QRHO) = qm(i+1,j,k,QRHO) + sourcr; - qm(i+1,j,k,QRHO) = amrex::max(qm(i+1,j,k,QRHO), lsmall_dens); - qm(i+1,j,k,QPRES) = qm(i+1,j,k,QPRES) + sourcp; - qm(i+1,j,k,QREINT) = qm(i+1,j,k,QREINT) + source; - } - - if (i >= vlo[0]) { - qp(i,j,k,QRHO) = qp(i,j,k,QRHO) + sourcr; - qp(i,j,k,QRHO) = amrex::max(qp(i,j,k,QRHO), lsmall_dens); - qp(i,j,k,QPRES) = qp(i,j,k,QPRES) + sourcp; - qp(i,j,k,QREINT) = qp(i,j,k,QREINT) + source; - } - } -#endif - }); } From f0647a6586b6766af1708ac7adedd4fd2b584c62 Mon Sep 17 00:00:00 2001 From: Michael Zingale Date: Tue, 30 May 2023 15:12:38 -0400 Subject: [PATCH 03/11] fix 3d compilation --- Source/hydro/trace_ppm.cpp | 6 ++++++ 1 file changed, 6 insertions(+) diff --git a/Source/hydro/trace_ppm.cpp b/Source/hydro/trace_ppm.cpp index 3f157305a9..9fa38b43c4 100644 --- a/Source/hydro/trace_ppm.cpp +++ b/Source/hydro/trace_ppm.cpp @@ -294,9 +294,11 @@ Castro::trace_ppm(const Box& bx, if (do_trace) { load_stencil(srcQ, idir, i, j, k, QRHO, s); +#if AMREX_SPACEDIM <= 2 if (idir == 0 && coord > 0) { add_geometric_rho_source(q_arr, dloga, i, j, k, s); } +#endif ppm_reconstruct(s, flat, sm, sp); ppm_int_profile(sm, sp, s[i0], un, cc, dtdx, Ip_src_rho, Im_src_rho); } @@ -338,9 +340,11 @@ Castro::trace_ppm(const Box& bx, if (do_trace) { load_stencil(srcQ, idir, i, j, k, QPRES, s); +#if AMREX_SPACEDIM <= 2 if (idir == 0 && coord > 0) { add_geometric_p_source(q_arr, qaux_arr, dloga, i, j, k, s); } +#endif ppm_reconstruct(s, flat, sm, sp); ppm_int_profile(sm, sp, s[i0], un, cc, dtdx, Ip_src_p, Im_src_p); } @@ -362,9 +366,11 @@ Castro::trace_ppm(const Box& bx, if (do_trace) { load_stencil(srcQ, idir, i, j, k, QREINT, s); +#if AMREX_SPACEDIM <= 2 if (idir == 0 && coord > 0) { add_geometric_rhoe_source(q_arr, dloga, i, j, k, s); } +#endif ppm_reconstruct(s, flat, sm, sp); ppm_int_profile(sm, sp, s[i0], un, cc, dtdx, Ip_src_rhoe, Im_src_rhoe); } From 8828c4bf1950ea37d5da36d4073af63eb2537bc0 Mon Sep 17 00:00:00 2001 From: Michael Zingale Date: Sat, 27 Jul 2024 10:46:45 -0400 Subject: [PATCH 04/11] fix pressure source --- Source/hydro/reconstruction.H | 52 +++++++++++++++++++++-------------- 1 file changed, 32 insertions(+), 20 deletions(-) diff --git a/Source/hydro/reconstruction.H b/Source/hydro/reconstruction.H index 65d6fc31dc..9a41710663 100644 --- a/Source/hydro/reconstruction.H +++ b/Source/hydro/reconstruction.H @@ -1,6 +1,8 @@ #ifndef CASTRO_RECONSTRUCTION_H #define CASTRO_RECONSTRUCTION_H +#include + namespace reconstruction { enum slope_indices { im2 = 0, @@ -14,9 +16,9 @@ namespace reconstruction { AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void -load_stencil(Array4 const& q_arr, const int idir, +load_stencil(amrex::Array4 const& q_arr, const int idir, const int i, const int j, const int k, const int ncomp, - Real* s) { + amrex::Real* s) { using namespace reconstruction; @@ -47,9 +49,11 @@ load_stencil(Array4 const& q_arr, const int idir, AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void -load_passive_stencil(Array4 const& U_arr, Array4 const& rho_inv_arr, const int idir, +load_passive_stencil(amrex::Array4 const& U_arr, + amrex::Array4 const& rho_inv_arr, + const int idir, const int i, const int j, const int k, const int ncomp, - Real* s) { + amrex::Real* s) { using namespace reconstruction; @@ -80,13 +84,16 @@ load_passive_stencil(Array4 const& U_arr, Array4 const& AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void -add_geometric_rho_source(Array4 const& q_arr, - Array4 const& dloga, +add_geometric_rho_source(amrex::Array4 const& q_arr, + amrex::Array4 const& dloga, const int i, const int j, const int k, - Real* s) { + amrex::Real* s) { using namespace reconstruction; + // this takes the form: -alpha rho u / r + // where alpha = 1 for cylindrical and 2 for spherical + // note: this is assumed to be working only in the x-direction s[im2] += -dloga(i-2,j,k) * q_arr(i-2,j,k,QRHO) * q_arr(i-2,j,k,QU); @@ -98,13 +105,16 @@ add_geometric_rho_source(Array4 const& q_arr, AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void -add_geometric_rhoe_source(Array4 const& q_arr, - Array4 const& dloga, +add_geometric_rhoe_source(amrex::Array4 const& q_arr, + amrex::Array4 const& dloga, const int i, const int j, const int k, - Real* s) { + amrex::Real* s) { using namespace reconstruction; + // this takes the form: -alpha (rho e + p) u / r + // where alpha = 1 for cylindrical and 2 for spherical + // note: this is assumed to be working only in the x-direction s[im2] += -dloga(i-2,j,k) * (q_arr(i-2,j,k,QREINT) + q_arr(i-2,j,k,QPRES)) * q_arr(i-2,j,k,QU); @@ -116,23 +126,25 @@ add_geometric_rhoe_source(Array4 const& q_arr, AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void -add_geometric_p_source(Array4 const& q_arr, - Array4 const& qaux_arr, - Array4 const& dloga, +add_geometric_p_source(amrex::Array4 const& q_arr, + amrex::Array4 const& qaux_arr, + amrex::Array4 const& dloga, const int i, const int j, const int k, - Real* s) { + amrex::Real* s) { using namespace reconstruction; + // this takes the form: -alpha Gamma1 p u / r + // where alpha = 1 for cylindrical and 2 for spherical + // note: this is assumed to be working only in the x-direction - s[im2] += -dloga(i-2,j,k) * q_arr(i-2,j,k,QRHO) * qaux_arr(i-2,j,k,QC) * q_arr(i-2,j,k,QU); - s[im1] += -dloga(i-1,j,k) * q_arr(i-1,j,k,QRHO) * qaux_arr(i-1,j,k,QC) * q_arr(i-1,j,k,QU); - s[i0] += -dloga(i,j,k) * q_arr(i,j,k,QRHO) * qaux_arr(i,j,k,QC) * q_arr(i,j,k,QU); - s[ip1] += -dloga(i+1,j,k) * q_arr(i+1,j,k,QRHO) * qaux_arr(i+1,j,k,QC) * q_arr(i+1,j,k,QU); - s[ip2] += -dloga(i+2,j,k) * q_arr(i+2,j,k,QRHO) * qaux_arr(i+2,j,k,QC) * q_arr(i+2,j,k,QU); + s[im2] += -dloga(i-2,j,k) * q_arr(i-2,j,k,QPRES) * qaux_arr(i-2,j,k,QGAMC) * q_arr(i-2,j,k,QU); + s[im1] += -dloga(i-1,j,k) * q_arr(i-1,j,k,QPRES) * qaux_arr(i-1,j,k,QGAMC) * q_arr(i-1,j,k,QU); + s[i0] += -dloga(i,j,k) * q_arr(i,j,k,QPRES) * qaux_arr(i,j,k,QGAMC) * q_arr(i,j,k,QU); + s[ip1] += -dloga(i+1,j,k) * q_arr(i+1,j,k,QPRES) * qaux_arr(i+1,j,k,QGAMC) * q_arr(i+1,j,k,QU); + s[ip2] += -dloga(i+2,j,k) * q_arr(i+2,j,k,QPRES) * qaux_arr(i+2,j,k,QGAMC) * q_arr(i+2,j,k,QU); } #endif - From 2c57bc80c13cbc00230ce90fce865a0ef1fb2107 Mon Sep 17 00:00:00 2001 From: Michael Zingale Date: Wed, 21 Aug 2024 10:23:23 -0400 Subject: [PATCH 05/11] remove unused var --- Source/hydro/trace_ppm.cpp | 9 --------- 1 file changed, 9 deletions(-) diff --git a/Source/hydro/trace_ppm.cpp b/Source/hydro/trace_ppm.cpp index 29faa9bd1e..0795500825 100644 --- a/Source/hydro/trace_ppm.cpp +++ b/Source/hydro/trace_ppm.cpp @@ -155,16 +155,9 @@ Castro::trace_ppm(const Box& bx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept { - Real cc = qaux_arr(i,j,k,QC); - -#if AMREX_SPACEDIM < 3 - Real csq = cc*cc; -#endif - Real un = q_arr(i,j,k,QUN); - // do the parabolic reconstruction and compute the // integrals under the characteristic waves Real s[nslp]; @@ -640,5 +633,3 @@ Castro::trace_ppm(const Box& bx, }); } - - From 8e914c38740c4675a5c7ff29a03d8e453dbaebbe Mon Sep 17 00:00:00 2001 From: Michael Zingale Date: Wed, 21 Aug 2024 11:37:03 -0400 Subject: [PATCH 06/11] update CI --- .../flame_wave/ci-benchmarks/grid_diag.out | 4 +- .../flame_wave/ci-benchmarks/species_diag.out | 4 +- .../ci-benchmarks/wdmerger_collision_2D.out | 46 +++++++++---------- 3 files changed, 27 insertions(+), 27 deletions(-) diff --git a/Exec/science/flame_wave/ci-benchmarks/grid_diag.out b/Exec/science/flame_wave/ci-benchmarks/grid_diag.out index 973d3f2aba..3bfc78ca85 100644 --- a/Exec/science/flame_wave/ci-benchmarks/grid_diag.out +++ b/Exec/science/flame_wave/ci-benchmarks/grid_diag.out @@ -1,5 +1,5 @@ # COLUMN 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 # TIMESTEP TIME MASS XMOM YMOM ZMOM ANG. MOM. X ANG. MOM. Y ANG. MOM. Z KIN. ENERGY INT. ENERGY GAS ENERGY GRAV. ENERGY TOTAL ENERGY CENTER OF MASS X-LOC CENTER OF MASS Y-LOC CENTER OF MASS Z-LOC CENTER OF MASS X-VEL CENTER OF MASS Y-VEL CENTER OF MASS Z-VEL MAXIMUM TEMPERATURE MAXIMUM DENSITY MAXIMUM T_S / T_E 0 0.0000000000000000e+00 1.4700685690736437e+20 0.0000000000000000e+00 0.0000000000000000e+00 0.0000000000000000e+00 0.0000000000000000e+00 0.0000000000000000e+00 0.0000000000000000e+00 0.0000000000000000e+00 2.9163021935541997e+37 2.9163021935541997e+37 0.0000000000000000e+00 2.9163021935541997e+37 2.7306641645125415e+04 6.8087525965685256e+02 0.0000000000000000e+00 0.0000000000000000e+00 0.0000000000000000e+00 0.0000000000000000e+00 1.3928678114307070e+09 3.0715027784567930e+07 0.0000000000000000e+00 - 1 5.9806047036494849e-08 1.4700700110495903e+20 3.3564221294760178e+23 7.2876252590129013e+23 1.2612378576902532e+20 3.1318392026636616e+23 -2.7141686393328650e+24 1.9066097490832445e+28 7.8605333893026802e+29 2.9163081863790971e+37 2.9163082649844082e+37 0.0000000000000000e+00 2.9163082649844082e+37 2.7306641717925890e+04 6.8087482283157362e+02 0.0000000000000000e+00 2.2831716205676648e+03 4.9573321027137572e+03 8.5794407627549896e-01 1.3929561430929687e+09 3.0715326643214259e+07 3.3873763068343888e-04 - 2 1.2559269877663918e-07 1.4700717750115456e+20 7.0484895627215180e+23 1.5304335700608418e+24 5.5619062084401319e+20 1.3810901029080806e+24 -1.1969287510263666e+25 4.0039912973182075e+28 3.0853296206485770e+30 2.9163145904134766e+37 2.9163148989464290e+37 0.0000000000000000e+00 2.9163148989464290e+37 2.7306641954244220e+04 6.8087461718998884e+02 0.0000000000000000e+00 4.7946567525018709e+03 1.0410604407725754e+04 3.7834249340624542e+00 1.3933944931271181e+09 3.0715562157661017e+07 3.3857586853705914e-04 + 1 5.9806047036494849e-08 1.4700700110495903e+20 3.3564221217849212e+23 7.2876252590226133e+23 1.2612378556916341e+20 3.1318391966955516e+23 -2.7141686347841742e+24 1.9066097493127676e+28 7.8605324800324878e+29 2.9163081863790706e+37 2.9163082649843738e+37 0.0000000000000000e+00 2.9163082649843738e+37 2.7306641717925955e+04 6.8087482283157317e+02 0.0000000000000000e+00 2.2831716153358752e+03 4.9573321027203638e+03 8.5794407491595881e-01 1.3929561431310942e+09 3.0715326643212341e+07 3.3873763041849706e-04 + 2 1.2559269877663918e-07 1.4700717750115462e+20 7.0484895445039834e+23 1.5304335700823335e+24 5.5619061962888426e+20 1.3810900989697918e+24 -1.1969287484132617e+25 4.0039912979259806e+28 3.0853275723280981e+30 2.9163145904135894e+37 2.9163148989463411e+37 0.0000000000000000e+00 2.9163148989463411e+37 2.7306641954244522e+04 6.8087461718998338e+02 0.0000000000000000e+00 4.7946567401095926e+03 1.0410604407871944e+04 3.7834249257966728e+00 1.3933944939675827e+09 3.0715562157563787e+07 3.3857586854903492e-04 diff --git a/Exec/science/flame_wave/ci-benchmarks/species_diag.out b/Exec/science/flame_wave/ci-benchmarks/species_diag.out index 6baaafc600..f3c2dcfe7a 100644 --- a/Exec/science/flame_wave/ci-benchmarks/species_diag.out +++ b/Exec/science/flame_wave/ci-benchmarks/species_diag.out @@ -1,5 +1,5 @@ # COLUMN 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 # TIMESTEP TIME Mass He4 Mass C12 Mass O16 Mass Ne20 Mass Mg24 Mass Si28 Mass S32 Mass Ar36 Mass Ca40 Mass Ti44 Mass Cr48 Mass Fe52 Mass Ni56 0 0.0000000000000000e+00 2.8142378112400895e-15 7.3932235330395157e-24 7.3932235330395157e-24 7.3932235330395157e-24 7.3932235330395157e-24 7.3932235330395157e-24 7.3932235330395157e-24 7.3932235330395157e-24 7.3932235330395157e-24 7.3932235330395157e-24 7.3932235330395157e-24 7.3932235330395157e-24 7.1117997526547418e-14 - 1 5.9806047036494849e-08 2.8142323856670296e-15 5.4329394513687282e-21 7.4194476167127106e-24 7.3934515289195991e-24 7.3935541950036108e-24 7.3933935990475538e-24 7.3932321468746301e-24 7.3932312019180565e-24 7.3932308328287436e-24 7.3932307869393029e-24 7.3932307850741563e-24 7.3932307849933102e-24 7.1118070045957091e-14 - 2 1.2559269877663918e-07 2.8142264165380905e-15 1.1401993461614438e-20 7.4929577845086585e-24 7.3944696393124482e-24 7.3939305113509452e-24 7.3935838899709814e-24 7.3932425308155010e-24 7.3932405326302037e-24 7.3932397568503172e-24 7.3932396603620439e-24 7.3932396564405536e-24 7.3932396562705448e-24 7.1118158758588142e-14 + 1 5.9806047036494849e-08 2.8142323856670446e-15 5.4329394398418374e-21 7.4194443455948811e-24 7.3934514829832856e-24 7.3935541949306949e-24 7.3933935989353088e-24 7.3932321468740865e-24 7.3932312019180345e-24 7.3932308328287333e-24 7.3932307869392985e-24 7.3932307850741519e-24 7.3932307849933044e-24 7.1118070045957103e-14 + 2 1.2559269877663918e-07 2.8142264165381050e-15 1.1401993450709496e-20 7.4929528148974785e-24 7.3944694290388374e-24 7.3939305042586412e-24 7.3935838894001731e-24 7.3932425308122493e-24 7.3932405326301655e-24 7.3932397568503098e-24 7.3932396603620498e-24 7.3932396564405566e-24 7.3932396562705492e-24 7.1118158758588142e-14 diff --git a/Exec/science/wdmerger/ci-benchmarks/wdmerger_collision_2D.out b/Exec/science/wdmerger/ci-benchmarks/wdmerger_collision_2D.out index 88fbaf4f88..50abd61af7 100644 --- a/Exec/science/wdmerger/ci-benchmarks/wdmerger_collision_2D.out +++ b/Exec/science/wdmerger/ci-benchmarks/wdmerger_collision_2D.out @@ -1,29 +1,29 @@ plotfile = plt00086 time = 1.25 variables minimum value maximum value - density 8.693611703e-05 19565534.299 - xmom -5.4964100651e+14 1.3559128302e+14 - ymom -2.5530096328e+15 2.5530122744e+15 + density 8.6936468335e-05 19568007.309 + xmom -5.4959005475e+14 1.3559838981e+14 + ymom -2.5540206053e+15 2.5540214496e+15 zmom 0 0 - rho_E 7.4982062146e+11 5.0669247219e+24 - rho_e 7.1077581849e+11 5.0640768326e+24 - Temp 242288.68588 1409652233.6 - rho_He4 8.693611703e-17 3.599903302 - rho_C12 3.4774446812e-05 7825956.6934 - rho_O16 5.2161670217e-05 11739149.75 - rho_Ne20 8.693611703e-17 181951.0571 - rho_Mg24 8.693611703e-17 1192.7969729 - rho_Si28 8.693611703e-17 6.6913702949 - rho_S32 8.693611703e-17 0.00019493291655 - rho_Ar36 8.693611703e-17 1.9565534609e-05 - rho_Ca40 8.693611703e-17 1.9565534331e-05 - rho_Ti44 8.693611703e-17 1.9565534308e-05 - rho_Cr48 8.693611703e-17 1.9565534308e-05 - rho_Fe52 8.693611703e-17 1.9565534308e-05 - rho_Ni56 8.693611703e-17 1.9565534308e-05 - phiGrav -5.8708033902e+17 -2.3375498341e+16 - grav_x -684991644 -51428.243166 - grav_y -739606241.84 739606820.44 + rho_E 7.4987846833e+11 5.0676543192e+24 + rho_e 7.1083104678e+11 5.0648015049e+24 + Temp 242292.44287 1409581841.1 + rho_He4 8.6936468335e-17 3.5973411848 + rho_C12 3.4774587333e-05 7826946.398 + rho_O16 5.2161881e-05 11740633.889 + rho_Ne20 8.6936468335e-17 181819.49413 + rho_Mg24 8.6936468335e-17 1190.7712386 + rho_Si28 8.6936468335e-17 6.6817951193 + rho_S32 8.6936468335e-17 0.00019451843176 + rho_Ar36 8.6936468335e-17 1.9568007618e-05 + rho_Ca40 8.6936468335e-17 1.9568007341e-05 + rho_Ti44 8.6936468335e-17 1.9568007318e-05 + rho_Cr48 8.6936468335e-17 1.9568007318e-05 + rho_Fe52 8.6936468335e-17 1.9568007318e-05 + rho_Ni56 8.6936468335e-17 1.9568007318e-05 + phiGrav -5.8709462562e+17 -2.3375498549e+16 + grav_x -685026429.13 -51428.265677 + grav_y -739654246.49 739654206.24 grav_z 0 0 - rho_enuc -2.7633982574e+12 7.6429034885e+23 + rho_enuc -4.7815621457e+12 7.6360058391e+23 From be168962f35c8b5c3b2d2d11a01433c2a11a1652 Mon Sep 17 00:00:00 2001 From: Zhi Date: Mon, 23 Sep 2024 16:44:27 -0400 Subject: [PATCH 07/11] update trace ppm to be compatible with spherical 2d --- Source/hydro/reconstruction.H | 68 ++++++++++++++++++++++++----------- Source/hydro/trace_ppm.cpp | 67 +++++++++++++++++++--------------- 2 files changed, 85 insertions(+), 50 deletions(-) diff --git a/Source/hydro/reconstruction.H b/Source/hydro/reconstruction.H index 9a41710663..50c2a95062 100644 --- a/Source/hydro/reconstruction.H +++ b/Source/hydro/reconstruction.H @@ -87,20 +87,27 @@ void add_geometric_rho_source(amrex::Array4 const& q_arr, amrex::Array4 const& dloga, const int i, const int j, const int k, - amrex::Real* s) { + const int ncomp , amrex::Real* s) { using namespace reconstruction; + // For idir == 0, i.e. r-direction: // this takes the form: -alpha rho u / r // where alpha = 1 for cylindrical and 2 for spherical + // note dloga(idir==0) = alpha/r - // note: this is assumed to be working only in the x-direction + // For idir == 1, + // i.e. theta-direction for spherical and z-direction for cylindrical + // this takes the form: -rho v cot(theta) / r for spherical and 0 for cylindrical + // note: dloga(idir==1) = cot(theta)/r for spherical and 0 for cylindrical. + + // ncomp should be QU for idir == 0 and QV for idir == 1 - s[im2] += -dloga(i-2,j,k) * q_arr(i-2,j,k,QRHO) * q_arr(i-2,j,k,QU); - s[im1] += -dloga(i-1,j,k) * q_arr(i-1,j,k,QRHO) * q_arr(i-1,j,k,QU); - s[i0] += -dloga(i,j,k) * q_arr(i,j,k,QRHO) * q_arr(i,j,k,QU); - s[ip1] += -dloga(i+1,j,k) * q_arr(i+1,j,k,QRHO) * q_arr(i+1,j,k,QU); - s[ip2] += -dloga(i+2,j,k) * q_arr(i+2,j,k,QRHO) * q_arr(i+2,j,k,QU); + s[im2] += -dloga(i-2,j,k) * q_arr(i-2,j,k,QRHO) * q_arr(i-2,j,k,ncomp); + s[im1] += -dloga(i-1,j,k) * q_arr(i-1,j,k,QRHO) * q_arr(i-1,j,k,ncomp); + s[i0] += -dloga(i,j,k) * q_arr(i,j,k,QRHO) * q_arr(i,j,k,ncomp); + s[ip1] += -dloga(i+1,j,k) * q_arr(i+1,j,k,QRHO) * q_arr(i+1,j,k,ncomp); + s[ip2] += -dloga(i+2,j,k) * q_arr(i+2,j,k,QRHO) * q_arr(i+2,j,k,ncomp); } AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE @@ -108,20 +115,27 @@ void add_geometric_rhoe_source(amrex::Array4 const& q_arr, amrex::Array4 const& dloga, const int i, const int j, const int k, - amrex::Real* s) { + const int ncomp, amrex::Real* s) { using namespace reconstruction; - // this takes the form: -alpha (rho e + p) u / r + // For idir == 0, i.e. r-direction: + // this takes the form: -alpha (rho e + p ) u / r // where alpha = 1 for cylindrical and 2 for spherical + // note dloga(idir==0) = alpha/r - // note: this is assumed to be working only in the x-direction + // For idir == 1, + // i.e. theta-direction for spherical and z-direction for cylindrical + // this takes the form: - (rho e + p) v cot(theta) / r for spherical and 0 for cylindrical + // note: dloga(idir==1) = cot(theta)/r for spherical and 0 for cylindrical. - s[im2] += -dloga(i-2,j,k) * (q_arr(i-2,j,k,QREINT) + q_arr(i-2,j,k,QPRES)) * q_arr(i-2,j,k,QU); - s[im1] += -dloga(i-1,j,k) * (q_arr(i-1,j,k,QREINT) + q_arr(i-1,j,k,QPRES)) * q_arr(i-1,j,k,QU); - s[i0] += -dloga(i,j,k) * (q_arr(i,j,k,QREINT) + q_arr(i,j,k,QPRES)) * q_arr(i,j,k,QU); - s[ip1] += -dloga(i+1,j,k) * (q_arr(i+1,j,k,QREINT) + q_arr(i+1,j,k,QPRES)) * q_arr(i+1,j,k,QU); - s[ip2] += -dloga(i+2,j,k) * (q_arr(i+2,j,k,QREINT) + q_arr(i+2,j,k,QPRES)) * q_arr(i+2,j,k,QU); + // ncomp should be QU for idir == 0 and QV for idir == 1 + + s[im2] += -dloga(i-2,j,k) * (q_arr(i-2,j,k,QREINT) + q_arr(i-2,j,k,QPRES)) * q_arr(i-2,j,k,ncomp); + s[im1] += -dloga(i-1,j,k) * (q_arr(i-1,j,k,QREINT) + q_arr(i-1,j,k,QPRES)) * q_arr(i-1,j,k,comp); + s[i0] += -dloga(i,j,k) * (q_arr(i,j,k,QREINT) + q_arr(i,j,k,QPRES)) * q_arr(i,j,k,ncomp); + s[ip1] += -dloga(i+1,j,k) * (q_arr(i+1,j,k,QREINT) + q_arr(i+1,j,k,QPRES)) * q_arr(i+1,j,k,ncomp); + s[ip2] += -dloga(i+2,j,k) * (q_arr(i+2,j,k,QREINT) + q_arr(i+2,j,k,QPRES)) * q_arr(i+2,j,k,ncomp); } AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE @@ -130,20 +144,32 @@ add_geometric_p_source(amrex::Array4 const& q_arr, amrex::Array4 const& qaux_arr, amrex::Array4 const& dloga, const int i, const int j, const int k, - amrex::Real* s) { + const int ncomp amrex::Real* s) { using namespace reconstruction; + // For idir == 0, i.e. r-direction: + // this takes the form: -alpha Gamma1 p u / r + // where alpha = 1 for cylindrical and 2 for spherical + // note: dloga(idir==0) = alpha/r + + // For idir == 1, + // i.e. theta-direction for spherical and z-direction for cylindrical + // this takes the form: - Gamma1 p v cot(theta) / r for spherical and 0 for cylindrical + // note: dloga(idir==1) = cot(theta)/r for spherical and 0 for cylindrical. + + // ncomp should be QU for idir == 0 and QV for idir == 1 + // this takes the form: -alpha Gamma1 p u / r // where alpha = 1 for cylindrical and 2 for spherical // note: this is assumed to be working only in the x-direction - s[im2] += -dloga(i-2,j,k) * q_arr(i-2,j,k,QPRES) * qaux_arr(i-2,j,k,QGAMC) * q_arr(i-2,j,k,QU); - s[im1] += -dloga(i-1,j,k) * q_arr(i-1,j,k,QPRES) * qaux_arr(i-1,j,k,QGAMC) * q_arr(i-1,j,k,QU); - s[i0] += -dloga(i,j,k) * q_arr(i,j,k,QPRES) * qaux_arr(i,j,k,QGAMC) * q_arr(i,j,k,QU); - s[ip1] += -dloga(i+1,j,k) * q_arr(i+1,j,k,QPRES) * qaux_arr(i+1,j,k,QGAMC) * q_arr(i+1,j,k,QU); - s[ip2] += -dloga(i+2,j,k) * q_arr(i+2,j,k,QPRES) * qaux_arr(i+2,j,k,QGAMC) * q_arr(i+2,j,k,QU); + s[im2] += -dloga(i-2,j,k) * q_arr(i-2,j,k,QPRES) * qaux_arr(i-2,j,k,QGAMC) * q_arr(i-2,j,k,ncomp); + s[im1] += -dloga(i-1,j,k) * q_arr(i-1,j,k,QPRES) * qaux_arr(i-1,j,k,QGAMC) * q_arr(i-1,j,k,ncomp); + s[i0] += -dloga(i,j,k) * q_arr(i,j,k,QPRES) * qaux_arr(i,j,k,QGAMC) * q_arr(i,j,k,ncomp); + s[ip1] += -dloga(i+1,j,k) * q_arr(i+1,j,k,QPRES) * qaux_arr(i+1,j,k,QGAMC) * q_arr(i+1,j,k,ncomp); + s[ip2] += -dloga(i+2,j,k) * q_arr(i+2,j,k,QPRES) * qaux_arr(i+2,j,k,QGAMC) * q_arr(i+2,j,k,ncomp); } diff --git a/Source/hydro/trace_ppm.cpp b/Source/hydro/trace_ppm.cpp index 0795500825..63b160e32e 100644 --- a/Source/hydro/trace_ppm.cpp +++ b/Source/hydro/trace_ppm.cpp @@ -55,10 +55,12 @@ Castro::trace_ppm(const Box& bx, const auto dx = geom.CellSizeArray(); + const auto problo = geom.ProbLoArray(); const int coord = geom.Coord(); Real hdt = 0.5_rt * dt; - Real dtdx = dt / dx[idir]; + Real dtdL = dt / dx[idir]; + Real dL = dx[idir]; #ifndef AMREX_USE_GPU auto lo = bx.loVect3d(); @@ -84,7 +86,8 @@ Castro::trace_ppm(const Box& bx, do_source_trace[n] = 0; // geometric source terms in r-direction need tracing - if (coord > 0 && idir == 0 && (n == QRHO || n == QPRES || n == QREINT)) { + if ((coord == 2 || (idir == 0 && coord == 1)) + && (n == QRHO || n == QPRES || n == QREINT)) { do_source_trace[n] = 1; continue; } @@ -154,6 +157,12 @@ Castro::trace_ppm(const Box& bx, amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept { + // Want dt/(rdtheta) instead of dt/dtheta for 2d Spherical + if (coord == 2 && idir == 1) { + Real r = problo[0] + static_cast(i + 0.5_rt) * dx[0]; + dL = r * dx[1]; + dtdL = dt / dL; + } Real cc = qaux_arr(i,j,k,QC); Real un = q_arr(i,j,k,QUN); @@ -196,7 +205,7 @@ Castro::trace_ppm(const Box& bx, load_stencil(q_arr, idir, i, j, k, QRHO, s); ppm_reconstruct(s, flat, sm, sp); - ppm_int_profile(sm, sp, s[i0], un, cc, dtdx, Ip_rho, Im_rho); + ppm_int_profile(sm, sp, s[i0], un, cc, dtdL, Ip_rho, Im_rho); // reconstruct normal velocity @@ -207,8 +216,8 @@ Castro::trace_ppm(const Box& bx, load_stencil(q_arr, idir, i, j, k, QUN, s); ppm_reconstruct(s, flat, sm, sp); - ppm_int_profile_single(sm, sp, s[i0], un-cc, dtdx, Ip_un_0, Im_un_0); - ppm_int_profile_single(sm, sp, s[i0], un+cc, dtdx, Ip_un_2, Im_un_2); + ppm_int_profile_single(sm, sp, s[i0], un-cc, dtdL, Ip_un_0, Im_un_0); + ppm_int_profile_single(sm, sp, s[i0], un+cc, dtdL, Ip_un_2, Im_un_2); // reconstruct pressure @@ -224,12 +233,12 @@ Castro::trace_ppm(const Box& bx, load_stencil(q_arr, idir, i, j, k, QRHO, trho); load_stencil(srcQ, idir, i, j, k, QUN, src); - ppm_reconstruct_pslope(trho, s, src, flat, dx[idir], sm, sp); + ppm_reconstruct_pslope(trho, s, src, flat, dL, sm, sp); } else { ppm_reconstruct(s, flat, sm, sp); } - ppm_int_profile(sm, sp, s[i0], un, cc, dtdx, Ip_p, Im_p); + ppm_int_profile(sm, sp, s[i0], un, cc, dtdL, Ip_p, Im_p); // reconstruct rho e @@ -238,7 +247,7 @@ Castro::trace_ppm(const Box& bx, load_stencil(q_arr, idir, i, j, k, QREINT, s); ppm_reconstruct(s, flat, sm, sp); - ppm_int_profile(sm, sp, s[i0], un, cc, dtdx, Ip_rhoe, Im_rhoe); + ppm_int_profile(sm, sp, s[i0], un, cc, dtdL, Ip_rhoe, Im_rhoe); // reconstruct transverse velocities @@ -249,11 +258,11 @@ Castro::trace_ppm(const Box& bx, load_stencil(q_arr, idir, i, j, k, QUT, s); ppm_reconstruct(s, flat, sm, sp); - ppm_int_profile_single(sm, sp, s[i0], un, dtdx, Ip_ut_1, Im_ut_1); + ppm_int_profile_single(sm, sp, s[i0], un, dtdL, Ip_ut_1, Im_ut_1); load_stencil(q_arr, idir, i, j, k, QUTT, s); ppm_reconstruct(s, flat, sm, sp); - ppm_int_profile_single(sm, sp, s[i0], un, dtdx, Ip_utt_1, Im_utt_1); + ppm_int_profile_single(sm, sp, s[i0], un, dtdL, Ip_utt_1, Im_utt_1); // gamma_c @@ -264,8 +273,8 @@ Castro::trace_ppm(const Box& bx, load_stencil(qaux_arr, idir, i, j, k, QGAMC, s); ppm_reconstruct(s, flat, sm, sp); - ppm_int_profile_single(sm, sp, s[i0], un-cc, dtdx, Ip_gc_0, Im_gc_0); - ppm_int_profile_single(sm, sp, s[i0], un+cc, dtdx, Ip_gc_2, Im_gc_2); + ppm_int_profile_single(sm, sp, s[i0], un-cc, dtdL, Ip_gc_0, Im_gc_0); + ppm_int_profile_single(sm, sp, s[i0], un+cc, dtdL, Ip_gc_2, Im_gc_2); // source terms -- we only trace if the terms in the stencil are non-zero @@ -279,7 +288,7 @@ Castro::trace_ppm(const Box& bx, #ifndef AMREX_USE_GPU do_trace = do_source_trace[QRHO]; #else - if (idir == 0 && coord > 0) { + if (coord == 2 || (idir == 0 && coord == 1)) { do_trace = 1; } else { do_trace = check_trace_source(srcQ, idir, i, j, k, QRHO); @@ -289,12 +298,12 @@ Castro::trace_ppm(const Box& bx, if (do_trace) { load_stencil(srcQ, idir, i, j, k, QRHO, s); #if AMREX_SPACEDIM <= 2 - if (idir == 0 && coord > 0) { - add_geometric_rho_source(q_arr, dloga, i, j, k, s); + if (coord == 2 || (idir == 0 && coord == 1)) { + add_geometric_rho_source(q_arr, dloga, i, j, k, QUN, s); } #endif ppm_reconstruct(s, flat, sm, sp); - ppm_int_profile(sm, sp, s[i0], un, cc, dtdx, Ip_src_rho, Im_src_rho); + ppm_int_profile(sm, sp, s[i0], un, cc, dtdL, Ip_src_rho, Im_src_rho); } // normal velocity @@ -313,8 +322,8 @@ Castro::trace_ppm(const Box& bx, if (do_trace) { load_stencil(srcQ, idir, i, j, k, QUN, s); ppm_reconstruct(s, flat, sm, sp); - ppm_int_profile_single(sm, sp, s[i0], un-cc, dtdx, Ip_src_un_0, Im_src_un_0); - ppm_int_profile_single(sm, sp, s[i0], un+cc, dtdx, Ip_src_un_2, Im_src_un_2); + ppm_int_profile_single(sm, sp, s[i0], un-cc, dtdL, Ip_src_un_0, Im_src_un_0); + ppm_int_profile_single(sm, sp, s[i0], un+cc, dtdL, Ip_src_un_2, Im_src_un_2); } // pressure @@ -325,7 +334,7 @@ Castro::trace_ppm(const Box& bx, #ifndef AMREX_USE_GPU do_trace = do_source_trace[QPRES]; #else - if (idir == 0 && coord > 0) { + if (coord == 2 || (idir == 0 && coord == 1)) { do_trace = 1; } else { do_trace = check_trace_source(srcQ, idir, i, j, k, QPRES); @@ -335,12 +344,12 @@ Castro::trace_ppm(const Box& bx, if (do_trace) { load_stencil(srcQ, idir, i, j, k, QPRES, s); #if AMREX_SPACEDIM <= 2 - if (idir == 0 && coord > 0) { - add_geometric_p_source(q_arr, qaux_arr, dloga, i, j, k, s); + if (coord == 2 || (idir == 0 && coord == 1)) { + add_geometric_p_source(q_arr, qaux_arr, dloga, i, j, k, QUN, s); } #endif ppm_reconstruct(s, flat, sm, sp); - ppm_int_profile(sm, sp, s[i0], un, cc, dtdx, Ip_src_p, Im_src_p); + ppm_int_profile(sm, sp, s[i0], un, cc, dtdL, Ip_src_p, Im_src_p); } // rho e @@ -351,7 +360,7 @@ Castro::trace_ppm(const Box& bx, #ifndef AMREX_USE_GPU do_trace = do_source_trace[QREINT]; #else - if (idir == 0 && coord > 0) { + if (coord == 2 || (idir == 0 && coord == 1)) { do_trace = 1; } else { do_trace = check_trace_source(srcQ, idir, i, j, k, QREINT); @@ -361,12 +370,12 @@ Castro::trace_ppm(const Box& bx, if (do_trace) { load_stencil(srcQ, idir, i, j, k, QREINT, s); #if AMREX_SPACEDIM <= 2 - if (idir == 0 && coord > 0) { - add_geometric_rhoe_source(q_arr, dloga, i, j, k, s); + if (coord == 2 || (idir == 0 && coord == 1)) { + add_geometric_rhoe_source(q_arr, dloga, i, j, k, QUN, s); } #endif ppm_reconstruct(s, flat, sm, sp); - ppm_int_profile(sm, sp, s[i0], un, cc, dtdx, Ip_src_rhoe, Im_src_rhoe); + ppm_int_profile(sm, sp, s[i0], un, cc, dtdL, Ip_src_rhoe, Im_src_rhoe); } // transverse velocities @@ -383,7 +392,7 @@ Castro::trace_ppm(const Box& bx, if (do_trace) { load_stencil(srcQ, idir, i, j, k, QUT, s); ppm_reconstruct(s, flat, sm, sp); - ppm_int_profile_single(sm, sp, s[i0], un, dtdx, Ip_src_ut_1, Im_src_ut_1); + ppm_int_profile_single(sm, sp, s[i0], un, dtdL, Ip_src_ut_1, Im_src_ut_1); } Real Ip_src_utt_1 = 0.0_rt; @@ -398,7 +407,7 @@ Castro::trace_ppm(const Box& bx, if (do_trace) { load_stencil(srcQ, idir, i, j, k, QUTT, s); ppm_reconstruct(s, flat, sm, sp); - ppm_int_profile_single(sm, sp, s[i0], un, dtdx, Ip_src_utt_1, Im_src_utt_1); + ppm_int_profile_single(sm, sp, s[i0], un, dtdL, Ip_src_utt_1, Im_src_utt_1); } @@ -416,7 +425,7 @@ Castro::trace_ppm(const Box& bx, load_passive_stencil(U_arr, rho_inv_arr, idir, i, j, k, nc, s); ppm_reconstruct(s, flat, sm, sp); - ppm_int_profile_single(sm, sp, s[i0], un, dtdx, Ip_passive, Im_passive); + ppm_int_profile_single(sm, sp, s[i0], un, dtdL, Ip_passive, Im_passive); // Plus state on face i From a2794fb7eb7208f5cfd6bce46de63e420301792e Mon Sep 17 00:00:00 2001 From: Zhi Date: Mon, 23 Sep 2024 16:49:20 -0400 Subject: [PATCH 08/11] missing comma --- Source/hydro/reconstruction.H | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/Source/hydro/reconstruction.H b/Source/hydro/reconstruction.H index 50c2a95062..858bdbce04 100644 --- a/Source/hydro/reconstruction.H +++ b/Source/hydro/reconstruction.H @@ -87,7 +87,7 @@ void add_geometric_rho_source(amrex::Array4 const& q_arr, amrex::Array4 const& dloga, const int i, const int j, const int k, - const int ncomp , amrex::Real* s) { + const int ncomp, amrex::Real* s) { using namespace reconstruction; @@ -144,7 +144,7 @@ add_geometric_p_source(amrex::Array4 const& q_arr, amrex::Array4 const& qaux_arr, amrex::Array4 const& dloga, const int i, const int j, const int k, - const int ncomp amrex::Real* s) { + const int ncomp, amrex::Real* s) { using namespace reconstruction; From 592fc9c576c3e16258b4c17b71caa99971d94420 Mon Sep 17 00:00:00 2001 From: Zhi Date: Mon, 23 Sep 2024 16:51:46 -0400 Subject: [PATCH 09/11] misspell --- Source/hydro/reconstruction.H | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Source/hydro/reconstruction.H b/Source/hydro/reconstruction.H index 858bdbce04..4315085532 100644 --- a/Source/hydro/reconstruction.H +++ b/Source/hydro/reconstruction.H @@ -132,7 +132,7 @@ add_geometric_rhoe_source(amrex::Array4 const& q_arr, // ncomp should be QU for idir == 0 and QV for idir == 1 s[im2] += -dloga(i-2,j,k) * (q_arr(i-2,j,k,QREINT) + q_arr(i-2,j,k,QPRES)) * q_arr(i-2,j,k,ncomp); - s[im1] += -dloga(i-1,j,k) * (q_arr(i-1,j,k,QREINT) + q_arr(i-1,j,k,QPRES)) * q_arr(i-1,j,k,comp); + s[im1] += -dloga(i-1,j,k) * (q_arr(i-1,j,k,QREINT) + q_arr(i-1,j,k,QPRES)) * q_arr(i-1,j,k,ncomp); s[i0] += -dloga(i,j,k) * (q_arr(i,j,k,QREINT) + q_arr(i,j,k,QPRES)) * q_arr(i,j,k,ncomp); s[ip1] += -dloga(i+1,j,k) * (q_arr(i+1,j,k,QREINT) + q_arr(i+1,j,k,QPRES)) * q_arr(i+1,j,k,ncomp); s[ip2] += -dloga(i+2,j,k) * (q_arr(i+2,j,k,QREINT) + q_arr(i+2,j,k,QPRES)) * q_arr(i+2,j,k,ncomp); From 12803d5cc52eb8faf60f6ec8416398e9edb59450 Mon Sep 17 00:00:00 2001 From: Zhi Date: Mon, 23 Sep 2024 16:56:42 -0400 Subject: [PATCH 10/11] fix compilation --- Source/hydro/trace_ppm.cpp | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/Source/hydro/trace_ppm.cpp b/Source/hydro/trace_ppm.cpp index 63b160e32e..6ae50ef777 100644 --- a/Source/hydro/trace_ppm.cpp +++ b/Source/hydro/trace_ppm.cpp @@ -59,8 +59,6 @@ Castro::trace_ppm(const Box& bx, const int coord = geom.Coord(); Real hdt = 0.5_rt * dt; - Real dtdL = dt / dx[idir]; - Real dL = dx[idir]; #ifndef AMREX_USE_GPU auto lo = bx.loVect3d(); @@ -157,6 +155,9 @@ Castro::trace_ppm(const Box& bx, amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept { + Real dtdL = dt / dx[idir]; + Real dL = dx[idir]; + // Want dt/(rdtheta) instead of dt/dtheta for 2d Spherical if (coord == 2 && idir == 1) { Real r = problo[0] + static_cast(i + 0.5_rt) * dx[0]; From 0ea7f1b7461293eaada403b991aba03c2b5e9dde Mon Sep 17 00:00:00 2001 From: Zhi Date: Mon, 23 Sep 2024 17:22:05 -0400 Subject: [PATCH 11/11] delete comment --- Source/hydro/reconstruction.H | 5 ----- 1 file changed, 5 deletions(-) diff --git a/Source/hydro/reconstruction.H b/Source/hydro/reconstruction.H index 4315085532..8f791a1ad5 100644 --- a/Source/hydro/reconstruction.H +++ b/Source/hydro/reconstruction.H @@ -160,11 +160,6 @@ add_geometric_p_source(amrex::Array4 const& q_arr, // ncomp should be QU for idir == 0 and QV for idir == 1 - // this takes the form: -alpha Gamma1 p u / r - // where alpha = 1 for cylindrical and 2 for spherical - - // note: this is assumed to be working only in the x-direction - s[im2] += -dloga(i-2,j,k) * q_arr(i-2,j,k,QPRES) * qaux_arr(i-2,j,k,QGAMC) * q_arr(i-2,j,k,ncomp); s[im1] += -dloga(i-1,j,k) * q_arr(i-1,j,k,QPRES) * qaux_arr(i-1,j,k,QGAMC) * q_arr(i-1,j,k,ncomp); s[i0] += -dloga(i,j,k) * q_arr(i,j,k,QPRES) * qaux_arr(i,j,k,QGAMC) * q_arr(i,j,k,ncomp);