From 61f0e38a990dec6449aad743b4603cd7962d935e Mon Sep 17 00:00:00 2001 From: Ann Almgren Date: Tue, 12 Mar 2024 09:31:55 -0700 Subject: [PATCH] fixes based on PR feedback --- Src/AmrCore/AMReX_InterpFaceRegister.H | 2 ++ Src/AmrCore/AMReX_InterpFaceRegister.cpp | 10 +++++----- Src/AmrCore/AMReX_Interp_C.H | 12 ++++++------ Src/AmrCore/AMReX_Interpolater.cpp | 1 - 4 files changed, 13 insertions(+), 12 deletions(-) diff --git a/Src/AmrCore/AMReX_InterpFaceRegister.H b/Src/AmrCore/AMReX_InterpFaceRegister.H index c54879bcaf6..42e38df29c0 100644 --- a/Src/AmrCore/AMReX_InterpFaceRegister.H +++ b/Src/AmrCore/AMReX_InterpFaceRegister.H @@ -74,6 +74,8 @@ private: Geometry m_crse_geom; + Box per_grown_domain; + Array m_fine_face_ba; Array m_crse_face_ba; Array m_face_mask; // crse/fine: 1, fine/fine & fine/physbc: 0 diff --git a/Src/AmrCore/AMReX_InterpFaceRegister.cpp b/Src/AmrCore/AMReX_InterpFaceRegister.cpp index 33de90c15b4..858b2c0292a 100644 --- a/Src/AmrCore/AMReX_InterpFaceRegister.cpp +++ b/Src/AmrCore/AMReX_InterpFaceRegister.cpp @@ -21,7 +21,7 @@ void InterpFaceRegister::define (BoxArray const& fba, DistributionMapping const& m_crse_geom = amrex::coarsen(m_fine_geom, m_ref_ratio); // We don't need to worry about face-based domain because this is only used in the tangential interpolation - Box per_grown_domain = m_crse_geom.Domain(); + per_grown_domain = m_crse_geom.Domain(); for (int dim = 0; dim < AMREX_SPACEDIM; dim++) { if (m_crse_geom.isPeriodic(dim)) { per_grown_domain.grow(dim,1); @@ -47,6 +47,7 @@ void InterpFaceRegister::define (BoxArray const& fba, DistributionMapping const& #ifdef AMREX_USE_GPU Vector > tags; #endif + #ifdef AMREX_USE_OMP #pragma omp parallel if (Gpu::notInLaunchRegion()) #endif @@ -99,7 +100,6 @@ namespace { Array4 slope; Array4 crse; Array4 mask; - Box per_grown_domain; AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Box box() const noexcept { return Box(mask); } }; @@ -159,7 +159,7 @@ InterpFaceRegister::interp (Array const& fine, // NOL { if (tag.mask(i,j,k)) { interp_face_reg(AMREX_D_DECL(i,j,k), rr, tag.fine, scomp, tag.crse, - tag.slope, ncomp, tag.per_grown_domain, idim); + tag.slope, ncomp, per_grown_domain, idim); } }); } else @@ -187,7 +187,7 @@ InterpFaceRegister::interp (Array const& fine, // NOL { if (mlo_arr(i,j,k)) { interp_face_reg(AMREX_D_DECL(i,j,k), rr, fine_arr, scomp, clo_arr, - slope_arr, ncomp, domlo, idim); + slope_arr, ncomp, per_grown_domain, idim); } }); slope.resize(hibx, ncomp); @@ -196,7 +196,7 @@ InterpFaceRegister::interp (Array const& fine, // NOL { if (mhi_arr(i,j,k)) { interp_face_reg(AMREX_D_DECL(i,j,k), rr, fine_arr, scomp, chi_arr, - slope_arr, ncomp, domhi, idim); + slope_arr, ncomp, per_grown_domain, idim); } }); } diff --git a/Src/AmrCore/AMReX_Interp_C.H b/Src/AmrCore/AMReX_Interp_C.H index 0b2a2a74bb0..01e920680c0 100644 --- a/Src/AmrCore/AMReX_Interp_C.H +++ b/Src/AmrCore/AMReX_Interp_C.H @@ -115,7 +115,7 @@ face_cons_linear_face_interp (int i, int j, int k, int n, Array4 const& fine, if (!mask || mask(ci, cj, ck, n)) { fine(i, j, k, n) = crse(ci, cj, ck, n); #if (AMREX_SPACEDIM >= 2) - if (cj != per_grown_domain.smallEnd(1) && cj != per_grown_domain.bigEnd(1) && ratio[1] > 1) { + if (cj > per_grown_domain.smallEnd(1) && cj < per_grown_domain.bigEnd(1) && ratio[1] > 1) { Real sfy = Real(1.0); Real dc = Real(0.5) * (crse(ci,cj+1,ck,n) - crse(ci,cj-1,ck,n)); Real df = Real(2.0) * (crse(ci,cj+1,ck,n) - crse(ci,cj ,ck,n)); @@ -131,7 +131,7 @@ face_cons_linear_face_interp (int i, int j, int k, int n, Array4 const& fine, fine(i,j,k,n) += yoff * slope * sfy; } // jc #if (AMREX_SPACEDIM == 3) - if (ck != per_grown_domain.smallEnd(2) && ck != per_grown_domain.bigEnd(2) && ratio[2] > 1) { + if (ck > per_grown_domain.smallEnd(2) && ck < per_grown_domain.bigEnd(2) && ratio[2] > 1) { Real sfz = Real(1.0); Real dc = Real(0.5) * (crse(ci,cj,ck+1,n) - crse(ci,cj,ck-1,n)); Real df = Real(2.0) * (crse(ci,cj,ck+1,n) - crse(ci,cj,ck ,n)); @@ -156,7 +156,7 @@ face_cons_linear_face_interp (int i, int j, int k, int n, Array4 const& fine, // Check solve mask to ensure we don't overwrite valid fine data. if (!mask || mask(ci, cj, ck, n)) { fine(i, j, k, n) = crse(ci, cj, ck, n); - if (ci != per_grown_domain.smallEnd(0) && ci != per_grown_domain.bigEnd(0) && ratio[0] > 1) { + if (ci > per_grown_domain.smallEnd(0) && ci < per_grown_domain.bigEnd(0) && ratio[0] > 1) { Real sfx = Real(1.0); Real dc = Real(0.5) * (crse(ci+1,cj,ck,n) - crse(ci-1,cj,ck,n)); Real df = Real(2.0) * (crse(ci+1,cj,ck,n) - crse(ci ,cj,ck,n)); @@ -172,7 +172,7 @@ face_cons_linear_face_interp (int i, int j, int k, int n, Array4 const& fine, fine(i,j,k,n) += xoff * slope * sfx; } // ci #if (AMREX_SPACEDIM == 3) - if (ck != per_grown_domain.smallEnd(2) && ck != per_grown_domain.bigEnd(2) && ratio[2] > 1) { + if (ck > per_grown_domain.smallEnd(2) && ck < per_grown_domain.bigEnd(2) && ratio[2] > 1) { Real sfz = Real(1.0); Real dc = Real(0.5) * (crse(ci,cj,ck+1,n) - crse(ci,cj,ck-1,n)); Real df = Real(2.0) * (crse(ci,cj,ck+1,n) - crse(ci,cj,ck ,n)); @@ -197,7 +197,7 @@ face_cons_linear_face_interp (int i, int j, int k, int n, Array4 const& fine, // Check solve mask to ensure we don't overwrite valid fine data. if (!mask || mask(ci, cj, ck, n)) { fine(i, j, k, n) = crse(ci, cj, ck, n); - if (ci != per_grown_domain.smallEnd(0) && ci != per_grown_domain.bigEnd(0) && ratio[0] > 1) { + if (ci > per_grown_domain.smallEnd(0) && ci < per_grown_domain.bigEnd(0) && ratio[0] > 1) { Real sfx = Real(1.0); Real dc = Real(0.5) * (crse(ci+1,cj,ck,n) - crse(ci-1,cj,ck,n)); Real df = Real(2.0) * (crse(ci+1,cj,ck,n) - crse(ci ,cj,ck,n)); @@ -212,7 +212,7 @@ face_cons_linear_face_interp (int i, int j, int k, int n, Array4 const& fine, Real xoff = (static_cast(i - ci*ratio[0]) + Real(0.5)) / Real(ratio[0]) - Real(0.5); fine(i,j,k,n) += xoff * slope * sfx; } // ci - if (cj != per_grown_domain.smallEnd(1) && cj != per_grown_domain.bigEnd(1) && ratio[1] > 1) { + if (cj > per_grown_domain.smallEnd(1) && cj < per_grown_domain.bigEnd(1) && ratio[1] > 1) { Real sfy = Real(1.0); Real dc = Real(0.5) * (crse(ci,cj+1,ck,n) - crse(ci,cj-1,ck,n)); Real df = Real(2.0) * (crse(ci,cj+1,ck,n) - crse(ci,cj ,ck,n)); diff --git a/Src/AmrCore/AMReX_Interpolater.cpp b/Src/AmrCore/AMReX_Interpolater.cpp index 532730c06b1..b5e855feb13 100644 --- a/Src/AmrCore/AMReX_Interpolater.cpp +++ b/Src/AmrCore/AMReX_Interpolater.cpp @@ -429,7 +429,6 @@ FaceConservativeLinear::CoarseBox (const Box& fine, const IntVect& ratio) for (int i = 0; i < AMREX_SPACEDIM; i++) { if (b.type(i) == IndexType::NODE) { - ng[i] = 0; if (b.type(i) == IndexType::NODE && b.length(i) < 2) { // Don't want degenerate boxes in nodal direction. b.growHi(i,1);