From e92cb5cdd62e65e502945cd8205d399ef274edfc Mon Sep 17 00:00:00 2001 From: Lorenzo Giacomel <47607756+lgiacome@users.noreply.github.com> Date: Tue, 27 Apr 2021 20:18:11 +0200 Subject: Staircased embedded boundaries in the YEE solver (#1881) * Added staircased embedded boundaris to the YEE solver * adding spherical resonating cavity test * adding functions for fields initialization * style adjustments * fixing tabs * fixed name of analysis script * fixed name of analysis script * fixed a few wrong preprocessor directives * workaround for missing boost * Revert "workaround for missing boost" This reverts commit 601f9eb2ec6f8c2100304379b2bea1c6cf9d1851. * another workaround for missing boost * getting rid of boost by depending on c++17 * Removed a few unused variables * adding USE_EB to addToCompileString for EB testing * removed tabs * fixing the inputs name for EB sphere test * shortened the test * zero padding the names of the images * adjusted two for loops * removed some unused variables * improving the fields initialization * removed the sphere test and implemented the cube test * fixed edges lengths computation and added comments * Fixed the case of all_regular geometries * fixing a bug that was breaking some tests * adding test folder * fixed the default values for the EB cube test * simplified the analysis script * fixed cubic resonator default results * inputting the plot file name from command line * fixing the diag name * Fixed a bug in edges initialization Co-authored-by: Remi Lehe * Adding comments to the staircased yee solver (thanks Remi) Co-authored-by: Remi Lehe * fixed the cube resonator test * removed an unused import Co-authored-by: Remi Lehe --- .../FieldSolver/FiniteDifferenceSolver/EvolveB.cpp | 26 +++++++++++++++++++--- 1 file changed, 23 insertions(+), 3 deletions(-) (limited to 'Source/FieldSolver/FiniteDifferenceSolver/EvolveB.cpp') diff --git a/Source/FieldSolver/FiniteDifferenceSolver/EvolveB.cpp b/Source/FieldSolver/FiniteDifferenceSolver/EvolveB.cpp index 2f66796f3..c24b869e0 100644 --- a/Source/FieldSolver/FiniteDifferenceSolver/EvolveB.cpp +++ b/Source/FieldSolver/FiniteDifferenceSolver/EvolveB.cpp @@ -25,6 +25,7 @@ using namespace amrex; void FiniteDifferenceSolver::EvolveB ( std::array< std::unique_ptr, 3 >& Bfield, std::array< std::unique_ptr, 3 > const& Efield, + std::array< std::unique_ptr, 3 > const& face_areas, int lev, amrex::Real const dt ) { // Select algorithm (The choice of algorithm is a runtime option, @@ -37,15 +38,15 @@ void FiniteDifferenceSolver::EvolveB ( #else if (m_do_nodal) { - EvolveBCartesian ( Bfield, Efield, lev, dt ); + EvolveBCartesian ( Bfield, Efield, face_areas, lev, dt ); } else if (m_fdtd_algo == MaxwellSolverAlgo::Yee) { - EvolveBCartesian ( Bfield, Efield, lev, dt ); + EvolveBCartesian ( Bfield, Efield, face_areas, lev, dt ); } else if (m_fdtd_algo == MaxwellSolverAlgo::CKC) { - EvolveBCartesian ( Bfield, Efield, lev, dt ); + EvolveBCartesian ( Bfield, Efield, face_areas, lev, dt ); #endif } else { @@ -61,6 +62,7 @@ template void FiniteDifferenceSolver::EvolveBCartesian ( std::array< std::unique_ptr, 3 >& Bfield, std::array< std::unique_ptr, 3 > const& Efield, + std::array< std::unique_ptr, 3 > const& face_areas, int lev, amrex::Real const dt ) { amrex::LayoutData* cost = WarpX::getCosts(lev); @@ -84,6 +86,12 @@ void FiniteDifferenceSolver::EvolveBCartesian ( Array4 const& Ey = Efield[1]->array(mfi); Array4 const& Ez = Efield[2]->array(mfi); +#ifdef AMREX_USE_EB + amrex::Array4 const& Sx = face_areas[0]->array(mfi); + amrex::Array4 const& Sy = face_areas[1]->array(mfi); + amrex::Array4 const& Sz = face_areas[2]->array(mfi); +#endif + // Extract stencil coefficients Real const * const AMREX_RESTRICT coefs_x = m_stencil_coefs_x.dataPtr(); int const n_coefs_x = m_stencil_coefs_x.size(); @@ -101,16 +109,28 @@ void FiniteDifferenceSolver::EvolveBCartesian ( amrex::ParallelFor(tbx, tby, tbz, [=] AMREX_GPU_DEVICE (int i, int j, int k){ +#ifdef AMREX_USE_EB + // Skip field push if this cell is fully covered by embedded boundaries + if (Sx(i, j, k) <= 0) return; +#endif Bx(i, j, k) += dt * T_Algo::UpwardDz(Ey, coefs_z, n_coefs_z, i, j, k) - dt * T_Algo::UpwardDy(Ez, coefs_y, n_coefs_y, i, j, k); }, [=] AMREX_GPU_DEVICE (int i, int j, int k){ +#ifdef AMREX_USE_EB + // Skip field push if this cell is fully covered by embedded boundaries + if (Sy(i, j, k) <= 0) return; +#endif By(i, j, k) += dt * T_Algo::UpwardDx(Ez, coefs_x, n_coefs_x, i, j, k) - dt * T_Algo::UpwardDz(Ex, coefs_z, n_coefs_z, i, j, k); }, [=] AMREX_GPU_DEVICE (int i, int j, int k){ +#ifdef AMREX_USE_EB + // Skip field push if this cell is fully covered by embedded boundaries + if (Sz(i, j, k) <= 0) return; +#endif Bz(i, j, k) += dt * T_Algo::UpwardDy(Ex, coefs_y, n_coefs_y, i, j, k) - dt * T_Algo::UpwardDx(Ey, coefs_x, n_coefs_x, i, j, k); } -- cgit v1.2.3