aboutsummaryrefslogtreecommitdiff
path: root/Source/FieldSolver/WarpXPushFieldsEM.cpp
diff options
context:
space:
mode:
Diffstat (limited to 'Source/FieldSolver/WarpXPushFieldsEM.cpp')
-rw-r--r--Source/FieldSolver/WarpXPushFieldsEM.cpp268
1 files changed, 35 insertions, 233 deletions
diff --git a/Source/FieldSolver/WarpXPushFieldsEM.cpp b/Source/FieldSolver/WarpXPushFieldsEM.cpp
index f63102b34..67fb40c01 100644
--- a/Source/FieldSolver/WarpXPushFieldsEM.cpp
+++ b/Source/FieldSolver/WarpXPushFieldsEM.cpp
@@ -140,81 +140,24 @@ void
WarpX::EvolveB (int lev, PatchType patch_type, amrex::Real a_dt)
{
+ // Evolve B field in regular cells
if (patch_type == PatchType::fine) {
m_fdtd_solver_fp[lev]->EvolveB( Bfield_fp[lev], Efield_fp[lev], a_dt );
} else {
m_fdtd_solver_cp[lev]->EvolveB( Bfield_cp[lev], Efield_cp[lev], a_dt );
}
- const int patch_level = (patch_type == PatchType::fine) ? lev : lev-1;
- const std::array<Real,3>& dx = WarpX::CellSize(patch_level);
- const Real dtsdx = a_dt/dx[0], dtsdy = a_dt/dx[1], dtsdz = a_dt/dx[2];
-
- if (do_pml && pml[lev]->ok())
- {
- const auto& pml_B = (patch_type == PatchType::fine) ? pml[lev]->GetB_fp() : pml[lev]->GetB_cp();
- const auto& pml_E = (patch_type == PatchType::fine) ? pml[lev]->GetE_fp() : pml[lev]->GetE_cp();
-
-#ifdef _OPENMP
-#pragma omp parallel if (Gpu::notInLaunchRegion())
-#endif
- for ( MFIter mfi(*pml_B[0], TilingIfNotGPU()); mfi.isValid(); ++mfi )
- {
- const Box& tbx = mfi.tilebox( pml_B[0]->ixType().toIntVect() );
- const Box& tby = mfi.tilebox( pml_B[1]->ixType().toIntVect() );
- const Box& tbz = mfi.tilebox( pml_B[2]->ixType().toIntVect() );
- auto const& pml_Bxfab = pml_B[0]->array(mfi);
- auto const& pml_Byfab = pml_B[1]->array(mfi);
- auto const& pml_Bzfab = pml_B[2]->array(mfi);
- auto const& pml_Exfab = pml_E[0]->array(mfi);
- auto const& pml_Eyfab = pml_E[1]->array(mfi);
- auto const& pml_Ezfab = pml_E[2]->array(mfi);
- if (WarpX::maxwell_fdtd_solver_id == 0) {
- amrex::ParallelFor(tbx, tby, tbz,
- [=] AMREX_GPU_DEVICE (int i, int j, int k) {
- warpx_push_pml_bx_yee(i,j,k,pml_Bxfab,pml_Eyfab,pml_Ezfab,
- dtsdy,dtsdz);
- },
- [=] AMREX_GPU_DEVICE (int i, int j, int k) {
- warpx_push_pml_by_yee(i,j,k,pml_Byfab,pml_Exfab,pml_Ezfab,
- dtsdx,dtsdz);
- },
- [=] AMREX_GPU_DEVICE (int i, int j, int k) {
- warpx_push_pml_bz_yee(i,j,k,pml_Bzfab,pml_Exfab,pml_Eyfab,
- dtsdx,dtsdy);
- });
- } else if (WarpX::maxwell_fdtd_solver_id == 1) {
- Real betaxy, betaxz, betayx, betayz, betazx, betazy;
- Real gammax, gammay, gammaz;
- Real alphax, alphay, alphaz;
- warpx_calculate_ckc_coefficients(dtsdx, dtsdy, dtsdz,
- betaxy, betaxz, betayx, betayz,
- betazx, betazy, gammax, gammay,
- gammaz, alphax, alphay, alphaz);
-
- amrex::ParallelFor(tbx, tby, tbz,
- [=] AMREX_GPU_DEVICE (int i, int j, int k) {
- warpx_push_pml_bx_ckc(i,j,k,pml_Bxfab,pml_Eyfab,pml_Ezfab,
- betaxy, betaxz, betayx, betayz,
- betazx, betazy, gammax, gammay,
- gammaz, alphax, alphay, alphaz);
- },
- [=] AMREX_GPU_DEVICE (int i, int j, int k) {
- warpx_push_pml_by_ckc(i,j,k,pml_Byfab,pml_Exfab,pml_Ezfab,
- betaxy, betaxz, betayx, betayz,
- betazx, betazy, gammax, gammay,
- gammaz, alphax, alphay, alphaz);
- },
- [=] AMREX_GPU_DEVICE (int i, int j, int k) {
- warpx_push_pml_bz_ckc(i,j,k,pml_Bzfab,pml_Exfab,pml_Eyfab,
- betaxy, betaxz, betayx, betayz,
- betazx, betazy, gammax, gammay,
- gammaz, alphax, alphay, alphaz);
- });
-
- }
+ // Evolve B field in PML cells
+ if (do_pml && pml[lev]->ok()) {
+ if (patch_type == PatchType::fine) {
+ m_fdtd_solver_fp[lev]->EvolveBPML(
+ pml[lev]->GetB_fp(), pml[lev]->GetE_fp(), a_dt );
+ } else {
+ m_fdtd_solver_cp[lev]->EvolveBPML(
+ pml[lev]->GetB_cp(), pml[lev]->GetE_cp(), a_dt );
}
}
+
}
void
@@ -240,7 +183,7 @@ WarpX::EvolveE (int lev, amrex::Real a_dt)
void
WarpX::EvolveE (int lev, PatchType patch_type, amrex::Real a_dt)
{
-
+ // Evolve E field in regular cells
if (patch_type == PatchType::fine) {
m_fdtd_solver_fp[lev]->EvolveE( Efield_fp[lev], Bfield_fp[lev],
current_fp[lev], F_fp[lev], a_dt );
@@ -249,142 +192,20 @@ WarpX::EvolveE (int lev, PatchType patch_type, amrex::Real a_dt)
current_cp[lev], F_cp[lev], a_dt );
}
- const Real mu_c2_dt = (PhysConst::mu0*PhysConst::c*PhysConst::c) * a_dt;
- const Real c2dt = (PhysConst::c*PhysConst::c) * a_dt;
-
- const int patch_level = (patch_type == PatchType::fine) ? lev : lev-1;
- const std::array<Real,3>& dx = WarpX::CellSize(patch_level);
- const Real dtsdx_c2 = c2dt/dx[0], dtsdy_c2 = c2dt/dx[1], dtsdz_c2 = c2dt/dx[2];
-
- MultiFab* F;
- if (patch_type == PatchType::fine)
- {
- F = F_fp[lev].get();
- }
- else if (patch_type == PatchType::coarse)
- {
- F = F_cp[lev].get();
- }
-
- if (do_pml && pml[lev]->ok())
- {
- const auto& pml_B = (patch_type == PatchType::fine) ? pml[lev]->GetB_fp() : pml[lev]->GetB_cp();
- const auto& pml_E = (patch_type == PatchType::fine) ? pml[lev]->GetE_fp() : pml[lev]->GetE_cp();
- const auto& pml_j = (patch_type == PatchType::fine) ? pml[lev]->Getj_fp() : pml[lev]->Getj_cp();
- const auto& pml_F = (patch_type == PatchType::fine) ? pml[lev]->GetF_fp() : pml[lev]->GetF_cp();
- const auto& sigba = (patch_type == PatchType::fine) ? pml[lev]->GetMultiSigmaBox_fp()
- : pml[lev]->GetMultiSigmaBox_cp();
-#ifdef _OPENMP
-#pragma omp parallel if (Gpu::notInLaunchRegion())
-#endif
- for ( MFIter mfi(*pml_E[0], TilingIfNotGPU()); mfi.isValid(); ++mfi )
- {
- const Box& tex = mfi.tilebox( pml_E[0]->ixType().toIntVect() );
- const Box& tey = mfi.tilebox( pml_E[1]->ixType().toIntVect() );
- const Box& tez = mfi.tilebox( pml_E[2]->ixType().toIntVect() );
-
- auto const& pml_Exfab = pml_E[0]->array(mfi);
- auto const& pml_Eyfab = pml_E[1]->array(mfi);
- auto const& pml_Ezfab = pml_E[2]->array(mfi);
- auto const& pml_Bxfab = pml_B[0]->array(mfi);
- auto const& pml_Byfab = pml_B[1]->array(mfi);
- auto const& pml_Bzfab = pml_B[2]->array(mfi);
-
- amrex::ParallelFor(tex, tey, tez,
- [=] AMREX_GPU_DEVICE (int i, int j, int k) {
- warpx_push_pml_ex_yee(i,j,k,pml_Exfab,pml_Byfab,pml_Bzfab,
- dtsdy_c2,dtsdz_c2);
- },
- [=] AMREX_GPU_DEVICE (int i, int j, int k) {
- warpx_push_pml_ey_yee(i,j,k,pml_Eyfab,pml_Bxfab,pml_Bzfab,
- dtsdx_c2,dtsdz_c2);
- },
- [=] AMREX_GPU_DEVICE (int i, int j, int k) {
- warpx_push_pml_ez_yee(i,j,k,pml_Ezfab,pml_Bxfab,pml_Byfab,
- dtsdx_c2,dtsdy_c2);
- });
-
- if (pml_has_particles) {
- // Update the E field in the PML, using the current
- // deposited by the particles in the PML
- auto const& pml_jxfab = pml_j[0]->array(mfi);
- auto const& pml_jyfab = pml_j[1]->array(mfi);
- auto const& pml_jzfab = pml_j[2]->array(mfi);
- const Real* sigmaj_x = sigba[mfi].sigma[0].data();
- const Real* sigmaj_y = sigba[mfi].sigma[1].data();
- const Real* sigmaj_z = sigba[mfi].sigma[2].data();
-
- int const x_lo = sigba[mfi].sigma[0].lo();
-#if (AMREX_SPACEDIM == 3)
- int const y_lo = sigba[mfi].sigma[1].lo();
- int const z_lo = sigba[mfi].sigma[2].lo();
-#else
- int const y_lo = 0;
- int const z_lo = sigba[mfi].sigma[1].lo();
-#endif
- amrex::ParallelFor( tex, tey, tez,
- [=] AMREX_GPU_DEVICE (int i, int j, int k) {
- push_ex_pml_current(i,j,k,
- pml_Exfab, pml_jxfab, sigmaj_y, sigmaj_z,
- y_lo, z_lo, mu_c2_dt);
- },
- [=] AMREX_GPU_DEVICE (int i, int j, int k) {
- push_ey_pml_current(i,j,k,
- pml_Eyfab, pml_jyfab, sigmaj_x, sigmaj_z,
- x_lo, z_lo, mu_c2_dt);
- },
- [=] AMREX_GPU_DEVICE (int i, int j, int k) {
- push_ez_pml_current(i,j,k,
- pml_Ezfab, pml_jzfab, sigmaj_x, sigmaj_y,
- x_lo, y_lo, mu_c2_dt);
- }
- );
- }
-
-
- if (pml_F)
- {
-
- auto const& pml_F_fab = pml_F->array(mfi);
-
- if (WarpX::maxwell_fdtd_solver_id == 0) {
-
- amrex::ParallelFor(tex, tey, tez,
- [=] AMREX_GPU_DEVICE (int i, int j, int k) {
- warpx_push_pml_ex_f_yee(i,j,k,pml_Exfab,pml_F_fab,dtsdx_c2);
- },
- [=] AMREX_GPU_DEVICE (int i, int j, int k) {
- warpx_push_pml_ey_f_yee(i,j,k,pml_Eyfab,pml_F_fab,dtsdy_c2);
- },
- [=] AMREX_GPU_DEVICE (int i, int j, int k) {
- warpx_push_pml_ez_f_yee(i,j,k,pml_Ezfab,pml_F_fab,dtsdz_c2);
- });
-
- } else if (WarpX::maxwell_fdtd_solver_id == 1) {
-
- Real betaxy, betaxz, betayx, betayz, betazx, betazy;
- Real gammax, gammay, gammaz;
- Real alphax, alphay, alphaz;
- warpx_calculate_ckc_coefficients(dtsdx_c2, dtsdy_c2, dtsdz_c2,
- betaxy, betaxz, betayx, betayz,
- betazx, betazy, gammax, gammay,
- gammaz, alphax, alphay, alphaz);
- amrex::ParallelFor(tex, tey, tez,
- [=] AMREX_GPU_DEVICE (int i, int j, int k) {
- warpx_push_pml_ex_f_ckc(i,j,k,pml_Exfab,pml_F_fab,
- alphax,betaxy,betaxz,gammax);
- },
- [=] AMREX_GPU_DEVICE (int i, int j, int k) {
- warpx_push_pml_ey_f_ckc(i,j,k,pml_Eyfab,pml_F_fab,
- alphay,betayx,betayz,gammay);
- },
- [=] AMREX_GPU_DEVICE (int i, int j, int k) {
- warpx_push_pml_ez_f_ckc(i,j,k,pml_Ezfab,pml_F_fab,
- alphaz,betazx,betazy,gammaz);
- });
-
- }
- }
+ // Evolve E field in PML cells
+ if (do_pml && pml[lev]->ok()) {
+ if (patch_type == PatchType::fine) {
+ m_fdtd_solver_fp[lev]->EvolveEPML(
+ pml[lev]->GetE_fp(), pml[lev]->GetB_fp(),
+ pml[lev]->Getj_fp(), pml[lev]->GetF_fp(),
+ pml[lev]->GetMultiSigmaBox_fp(),
+ a_dt, pml_has_particles );
+ } else {
+ m_fdtd_solver_cp[lev]->EvolveEPML(
+ pml[lev]->GetE_cp(), pml[lev]->GetB_cp(),
+ pml[lev]->Getj_cp(), pml[lev]->GetF_cp(),
+ pml[lev]->GetMultiSigmaBox_cp(),
+ a_dt, pml_has_particles );
}
}
}
@@ -418,6 +239,7 @@ WarpX::EvolveF (int lev, PatchType patch_type, amrex::Real a_dt, DtType a_dt_typ
const int rhocomp = (a_dt_type == DtType::FirstHalf) ? 0 : 1;
+ // Evolve F field in regular cells
if (patch_type == PatchType::fine) {
m_fdtd_solver_fp[lev]->EvolveF( F_fp[lev], Efield_fp[lev],
rho_fp[lev], rhocomp, a_dt );
@@ -426,37 +248,17 @@ WarpX::EvolveF (int lev, PatchType patch_type, amrex::Real a_dt, DtType a_dt_typ
rho_cp[lev], rhocomp, a_dt );
}
- const int patch_level = (patch_type == PatchType::fine) ? lev : lev-1;
- const auto& dx = WarpX::CellSize(patch_level);
- const std::array<Real,3> dtsdx {a_dt/dx[0], a_dt/dx[1], a_dt/dx[2]};
-
- if (do_pml && pml[lev]->ok())
- {
- const auto& pml_F = (patch_type == PatchType::fine) ? pml[lev]->GetF_fp() : pml[lev]->GetF_cp();
- const auto& pml_E = (patch_type == PatchType::fine) ? pml[lev]->GetE_fp() : pml[lev]->GetE_cp();
-
-#ifdef _OPENMP
-#pragma omp parallel if (Gpu::notInLaunchRegion())
-#endif
- for ( MFIter mfi(*pml_F, TilingIfNotGPU()); mfi.isValid(); ++mfi )
- {
- const Box& bx = mfi.tilebox();
-
- auto const& pml_F_fab = pml_F->array(mfi);
- auto const& pml_Exfab = pml_E[0]->array(mfi);
- auto const& pml_Eyfab = pml_E[1]->array(mfi);
- auto const& pml_Ezfab = pml_E[2]->array(mfi);
-
- amrex::ParallelFor(bx,
- [=] AMREX_GPU_DEVICE (int i, int j, int k)
- {
- warpx_push_pml_F(i, j, k, pml_F_fab, pml_Exfab,
- pml_Eyfab, pml_Ezfab,
- dtsdx[0], dtsdx[1], dtsdx[2]);
- });
-
+ // Evolve F field in PML cells
+ if (do_pml && pml[lev]->ok()) {
+ if (patch_type == PatchType::fine) {
+ m_fdtd_solver_fp[lev]->EvolveFPML(
+ pml[lev]->GetF_fp(), pml[lev]->GetE_fp(), a_dt );
+ } else {
+ m_fdtd_solver_cp[lev]->EvolveFPML(
+ pml[lev]->GetF_cp(), pml[lev]->GetE_cp(), a_dt );
}
}
+
}
#ifdef WARPX_DIM_RZ
generated by cgit v1.2.3 (git 2.46.0) at 2025-08-11 12:55:47 +0000