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/EvolveE.cpp | 26 +++++++++++++++++++--- 1 file changed, 23 insertions(+), 3 deletions(-) (limited to 'Source/FieldSolver/FiniteDifferenceSolver/EvolveE.cpp') diff --git a/Source/FieldSolver/FiniteDifferenceSolver/EvolveE.cpp b/Source/FieldSolver/FiniteDifferenceSolver/EvolveE.cpp index cd1e3bdf9..3cb433e97 100644 --- a/Source/FieldSolver/FiniteDifferenceSolver/EvolveE.cpp +++ b/Source/FieldSolver/FiniteDifferenceSolver/EvolveE.cpp @@ -28,6 +28,7 @@ void FiniteDifferenceSolver::EvolveE ( std::array< std::unique_ptr, 3 >& Efield, std::array< std::unique_ptr, 3 > const& Bfield, std::array< std::unique_ptr, 3 > const& Jfield, + std::array< std::unique_ptr, 3 > const& edge_lengths, std::unique_ptr const& Ffield, int lev, amrex::Real const dt ) { @@ -41,15 +42,15 @@ void FiniteDifferenceSolver::EvolveE ( #else if (m_do_nodal) { - EvolveECartesian ( Efield, Bfield, Jfield, Ffield, lev, dt ); + EvolveECartesian ( Efield, Bfield, Jfield, edge_lengths, Ffield, lev, dt ); } else if (m_fdtd_algo == MaxwellSolverAlgo::Yee) { - EvolveECartesian ( Efield, Bfield, Jfield, Ffield, lev, dt ); + EvolveECartesian ( Efield, Bfield, Jfield, edge_lengths, Ffield, lev, dt ); } else if (m_fdtd_algo == MaxwellSolverAlgo::CKC) { - EvolveECartesian ( Efield, Bfield, Jfield, Ffield, lev, dt ); + EvolveECartesian ( Efield, Bfield, Jfield, edge_lengths, Ffield, lev, dt ); #endif } else { @@ -66,6 +67,7 @@ void FiniteDifferenceSolver::EvolveECartesian ( std::array< std::unique_ptr, 3 >& Efield, std::array< std::unique_ptr, 3 > const& Bfield, std::array< std::unique_ptr, 3 > const& Jfield, + std::array< std::unique_ptr, 3 > const& edge_lengths, std::unique_ptr const& Ffield, int lev, amrex::Real const dt ) { @@ -94,6 +96,12 @@ void FiniteDifferenceSolver::EvolveECartesian ( Array4 const& jy = Jfield[1]->array(mfi); Array4 const& jz = Jfield[2]->array(mfi); +#ifdef AMREX_USE_EB + amrex::Array4 const& lx = edge_lengths[0]->array(mfi); + amrex::Array4 const& ly = edge_lengths[1]->array(mfi); + amrex::Array4 const& lz = edge_lengths[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(); @@ -111,6 +119,10 @@ void FiniteDifferenceSolver::EvolveECartesian ( amrex::ParallelFor(tex, tey, tez, [=] 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 (lx(i, j, k) <= 0) return; +#endif Ex(i, j, k) += c2 * dt * ( - T_Algo::DownwardDz(By, coefs_z, n_coefs_z, i, j, k) + T_Algo::DownwardDy(Bz, coefs_y, n_coefs_y, i, j, k) @@ -118,6 +130,10 @@ void FiniteDifferenceSolver::EvolveECartesian ( }, [=] 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 (ly(i,j,k) <= 0) return; +#endif Ey(i, j, k) += c2 * dt * ( - T_Algo::DownwardDx(Bz, coefs_x, n_coefs_x, i, j, k) + T_Algo::DownwardDz(Bx, coefs_z, n_coefs_z, i, j, k) @@ -125,6 +141,10 @@ void FiniteDifferenceSolver::EvolveECartesian ( }, [=] 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 (lz(i,j,k) <= 0) return; +#endif Ez(i, j, k) += c2 * dt * ( - T_Algo::DownwardDy(Bx, coefs_y, n_coefs_y, i, j, k) + T_Algo::DownwardDx(By, coefs_x, n_coefs_x, i, j, k) -- cgit v1.2.3