From c69297d58a75735d7995c27c1ba95a708fe9bc07 Mon Sep 17 00:00:00 2001 From: Edoardo Zoni <59625522+EZoni@users.noreply.github.com> Date: Wed, 9 Nov 2022 12:36:41 -0800 Subject: 2D/RZ Embedded Boundaries Bug Fix (#3510) * Run embedded_boundary_cube_2d in debug mode * Fix bug * Use `lbound`/`ubound` instead of `begin`/`end` * Run embedded_boundary_cube_2d in release mode --- Source/Initialization/WarpXInitData.cpp | 13 ++++++++++++- 1 file changed, 12 insertions(+), 1 deletion(-) (limited to 'Source/Initialization/WarpXInitData.cpp') diff --git a/Source/Initialization/WarpXInitData.cpp b/Source/Initialization/WarpXInitData.cpp index 3a5e9ed24..ab421f106 100644 --- a/Source/Initialization/WarpXInitData.cpp +++ b/Source/Initialization/WarpXInitData.cpp @@ -933,6 +933,7 @@ WarpX::InitializeExternalFieldsOnGridUsingParser ( amrex::IntVect x_nodal_flag = mfx->ixType().toIntVect(); amrex::IntVect y_nodal_flag = mfy->ixType().toIntVect(); amrex::IntVect z_nodal_flag = mfz->ixType().toIntVect(); + for ( MFIter mfi(*mfx, TilingIfNotGPU()); mfi.isValid(); ++mfi) { const amrex::Box& tbx = mfi.tilebox( x_nodal_flag, mfx->nGrowVect() ); @@ -951,6 +952,13 @@ WarpX::InitializeExternalFieldsOnGridUsingParser ( amrex::Array4 const& Sy = face_areas[1]->array(mfi); amrex::Array4 const& Sz = face_areas[2]->array(mfi); +#if defined(WARPX_DIM_XZ) || defined(WARPX_DIM_RZ) + const amrex::Dim3 lx_lo = amrex::lbound(lx); + const amrex::Dim3 lx_hi = amrex::ubound(lx); + const amrex::Dim3 lz_lo = amrex::lbound(lz); + const amrex::Dim3 lz_hi = amrex::ubound(lz); +#endif + #if defined(WARPX_DIM_XZ) || defined(WARPX_DIM_RZ) amrex::ignore_unused(ly, Sx, Sz); #elif defined(WARPX_DIM_1D_Z) @@ -1001,7 +1009,10 @@ WarpX::InitializeExternalFieldsOnGridUsingParser ( if((field=='E' and ly(i, j, k)<=0) or (field=='B' and Sy(i, j, k)<=0)) return; #elif defined(WARPX_DIM_XZ) || defined(WARPX_DIM_RZ) //In XZ and RZ Ey is associated with a mesh node, so we need to check if the mesh node is covered - if((field=='E' and (lx(i, j, k)<=0 || lx(i-1, j, k)<=0 || lz(i, j, k)<=0 || lz(i, j-1, k)<=0)) or + if((field=='E' and (lx(std::min(i , lx_hi.x), std::min(j , lx_hi.y), k)<=0 + || lx(std::max(i-1, lx_lo.x), std::min(j , lx_hi.y), k)<=0 + || lz(std::min(i , lz_hi.x), std::min(j , lz_hi.y), k)<=0 + || lz(std::min(i , lz_hi.x), std::max(j-1, lz_lo.y), k)<=0)) or (field=='B' and Sy(i,j,k)<=0)) return; #endif #endif -- cgit v1.2.3