aboutsummaryrefslogtreecommitdiff
path: root/Source/WarpXEvolve.cpp
diff options
context:
space:
mode:
Diffstat (limited to 'Source/WarpXEvolve.cpp')
-rw-r--r--Source/WarpXEvolve.cpp152
1 files changed, 76 insertions, 76 deletions
diff --git a/Source/WarpXEvolve.cpp b/Source/WarpXEvolve.cpp
index d142e175b..46dd2f25a 100644
--- a/Source/WarpXEvolve.cpp
+++ b/Source/WarpXEvolve.cpp
@@ -63,16 +63,16 @@ WarpX::EvolveEM (int numsteps)
if (step > 0 && (step+1) % load_balance_int == 0)
{
LoadBalance();
- // Reset the costs to 0
- for (int lev = 0; lev <= finest_level; ++lev) {
- costs[lev]->setVal(0.0);
- }
+ // Reset the costs to 0
+ for (int lev = 0; lev <= finest_level; ++lev) {
+ costs[lev]->setVal(0.0);
+ }
}
for (int lev = 0; lev <= finest_level; ++lev) {
- // Perform running average of the costs
- // (Giving more importance to most recent costs)
- (*costs[lev].get()).mult( (1. - 2./load_balance_int) );
+ // Perform running average of the costs
+ // (Giving more importance to most recent costs)
+ (*costs[lev].get()).mult( (1. - 2./load_balance_int) );
}
}
@@ -89,12 +89,12 @@ WarpX::EvolveEM (int numsteps)
*Efield_aux[lev][0],*Efield_aux[lev][1],*Efield_aux[lev][2],
*Bfield_aux[lev][0],*Bfield_aux[lev][1],*Bfield_aux[lev][2]);
}
- is_synchronized = false;
+ is_synchronized = false;
} else {
// Beyond one step, we have E^{n} and B^{n}.
// Particles have p^{n-1/2} and x^{n}.
UpdateAuxilaryData();
- }
+ }
// Push particle from x^{n} to x^{n+1}
// from p^{n-1/2} to p^{n+1/2}
@@ -146,7 +146,7 @@ WarpX::EvolveEM (int numsteps)
mypc->PushP(lev, 0.5*dt[0],
*Efield_aux[lev][0],*Efield_aux[lev][1],*Efield_aux[lev][2],
*Bfield_aux[lev][0],*Bfield_aux[lev][1],*Bfield_aux[lev][2]);
- }
+ }
is_synchronized = true;
}
#ifdef WARPX_USE_PY
@@ -157,14 +157,14 @@ WarpX::EvolveEM (int numsteps)
++istep[lev];
}
- cur_time += dt[0];
+ cur_time += dt[0];
bool to_make_plot = (plot_int > 0) && ((step+1) % plot_int == 0);
bool move_j = is_synchronized || to_make_plot;
// If is_synchronized we need to shift j too so that next step we can evolve E by dt/2.
// We might need to move j because we are going to make a plotfile.
- MoveWindow(move_j);
+ MoveWindow(move_j);
if (max_level == 0) {
mypc->RedistributeLocal();
@@ -181,10 +181,10 @@ WarpX::EvolveEM (int numsteps)
<< " s; This step = " << walltime_end_step-walltime_beg_step
<< " s; Avg. per step = " << walltime/(step+1) << " s\n";
- // sync up time
- for (int i = 0; i <= max_level; ++i) {
- t_new[i] = cur_time;
- }
+ // sync up time
+ for (int i = 0; i <= max_level; ++i) {
+ t_new[i] = cur_time;
+ }
if (do_boosted_frame_diagnostic) {
std::unique_ptr<MultiFab> cell_centered_data = nullptr;
@@ -194,7 +194,7 @@ WarpX::EvolveEM (int numsteps)
myBFD->writeLabFrameData(cell_centered_data.get(), *mypc, geom[0], cur_time, dt[0]);
}
- if (to_make_plot)
+ if (to_make_plot)
{
FillBoundaryE();
FillBoundaryB();
@@ -206,24 +206,24 @@ WarpX::EvolveEM (int numsteps)
*Bfield_aux[lev][0],*Bfield_aux[lev][1],*Bfield_aux[lev][2]);
}
- last_plot_file_step = step+1;
- WritePlotFile();
- }
+ last_plot_file_step = step+1;
+ WritePlotFile();
+ }
- if (check_int > 0 && (step+1) % check_int == 0) {
- last_check_file_step = step+1;
- WriteCheckPointFile();
- }
+ if (check_int > 0 && (step+1) % check_int == 0) {
+ last_check_file_step = step+1;
+ WriteCheckPointFile();
+ }
- if (cur_time >= stop_time - 1.e-3*dt[0]) {
- max_time_reached = true;
- break;
- }
+ if (cur_time >= stop_time - 1.e-3*dt[0]) {
+ max_time_reached = true;
+ break;
+ }
#ifdef WARPX_USE_PY
if (warpx_py_afterstep) warpx_py_afterstep();
#endif
- // End loop on time steps
+ // End loop on time steps
}
if (plot_int > 0 && istep[0] > last_plot_file_step && (max_time_reached || istep[0] >= max_step))
@@ -238,11 +238,11 @@ WarpX::EvolveEM (int numsteps)
*Bfield_aux[lev][0],*Bfield_aux[lev][1],*Bfield_aux[lev][2]);
}
- WritePlotFile();
+ WritePlotFile();
}
if (check_int > 0 && istep[0] > last_check_file_step && (max_time_reached || istep[0] >= max_step)) {
- WriteCheckPointFile();
+ WriteCheckPointFile();
}
if (do_boosted_frame_diagnostic) {
@@ -605,65 +605,65 @@ WarpX::EvolveF (int lev, PatchType patch_type, Real dt, DtType dt_type)
void
WarpX::DampPML ()
{
- for (int lev = 0; lev <= finest_level; ++lev) {
- DampPML(lev);
- }
+ for (int lev = 0; lev <= finest_level; ++lev) {
+ DampPML(lev);
+ }
}
void
WarpX::DampPML (int lev)
{
- DampPML(lev, PatchType::fine);
- if (lev > 0) DampPML(lev, PatchType::coarse);
+ DampPML(lev, PatchType::fine);
+ if (lev > 0) DampPML(lev, PatchType::coarse);
}
void
WarpX::DampPML (int lev, PatchType patch_type)
{
- if (!do_pml) return;
+ if (!do_pml) return;
- BL_PROFILE("WarpX::DampPML()");
+ BL_PROFILE("WarpX::DampPML()");
- if (pml[lev]->ok())
+ if (pml[lev]->ok())
{
- const auto& pml_E = (patch_type == PatchType::fine) ? pml[lev]->GetE_fp() : pml[lev]->GetE_cp();
- const auto& pml_B = (patch_type == PatchType::fine) ? pml[lev]->GetB_fp() : pml[lev]->GetB_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();
+ const auto& pml_E = (patch_type == PatchType::fine) ? pml[lev]->GetE_fp() : pml[lev]->GetE_cp();
+ const auto& pml_B = (patch_type == PatchType::fine) ? pml[lev]->GetB_fp() : pml[lev]->GetB_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
#endif
- for ( MFIter mfi(*pml_E[0],true); mfi.isValid(); ++mfi )
+ for ( MFIter mfi(*pml_E[0],true); mfi.isValid(); ++mfi )
{
- const Box& tex = mfi.tilebox(Ex_nodal_flag);
- const Box& tey = mfi.tilebox(Ey_nodal_flag);
- const Box& tez = mfi.tilebox(Ez_nodal_flag);
- const Box& tbx = mfi.tilebox(Bx_nodal_flag);
- const Box& tby = mfi.tilebox(By_nodal_flag);
- const Box& tbz = mfi.tilebox(Bz_nodal_flag);
-
- WRPX_DAMP_PML(tex.loVect(), tex.hiVect(),
- tey.loVect(), tey.hiVect(),
- tez.loVect(), tez.hiVect(),
- tbx.loVect(), tbx.hiVect(),
- tby.loVect(), tby.hiVect(),
- tbz.loVect(), tbz.hiVect(),
- BL_TO_FORTRAN_3D((*pml_E[0])[mfi]),
- BL_TO_FORTRAN_3D((*pml_E[1])[mfi]),
- BL_TO_FORTRAN_3D((*pml_E[2])[mfi]),
- BL_TO_FORTRAN_3D((*pml_B[0])[mfi]),
- BL_TO_FORTRAN_3D((*pml_B[1])[mfi]),
- BL_TO_FORTRAN_3D((*pml_B[2])[mfi]),
- WRPX_PML_TO_FORTRAN(sigba[mfi]));
+ const Box& tex = mfi.tilebox(Ex_nodal_flag);
+ const Box& tey = mfi.tilebox(Ey_nodal_flag);
+ const Box& tez = mfi.tilebox(Ez_nodal_flag);
+ const Box& tbx = mfi.tilebox(Bx_nodal_flag);
+ const Box& tby = mfi.tilebox(By_nodal_flag);
+ const Box& tbz = mfi.tilebox(Bz_nodal_flag);
- if (pml_F) {
- const Box& tnd = mfi.nodaltilebox();
- WRPX_DAMP_PML_F(tnd.loVect(), tnd.hiVect(),
- BL_TO_FORTRAN_3D((*pml_F)[mfi]),
+ WRPX_DAMP_PML(tex.loVect(), tex.hiVect(),
+ tey.loVect(), tey.hiVect(),
+ tez.loVect(), tez.hiVect(),
+ tbx.loVect(), tbx.hiVect(),
+ tby.loVect(), tby.hiVect(),
+ tbz.loVect(), tbz.hiVect(),
+ BL_TO_FORTRAN_3D((*pml_E[0])[mfi]),
+ BL_TO_FORTRAN_3D((*pml_E[1])[mfi]),
+ BL_TO_FORTRAN_3D((*pml_E[2])[mfi]),
+ BL_TO_FORTRAN_3D((*pml_B[0])[mfi]),
+ BL_TO_FORTRAN_3D((*pml_B[1])[mfi]),
+ BL_TO_FORTRAN_3D((*pml_B[2])[mfi]),
WRPX_PML_TO_FORTRAN(sigba[mfi]));
- }
+
+ if (pml_F) {
+ const Box& tnd = mfi.nodaltilebox();
+ WRPX_DAMP_PML_F(tnd.loVect(), tnd.hiVect(),
+ BL_TO_FORTRAN_3D((*pml_F)[mfi]),
+ WRPX_PML_TO_FORTRAN(sigba[mfi]));
+ }
}
}
}
@@ -699,18 +699,18 @@ WarpX::ComputeDt ()
Real deltat = 0.;
if (maxwell_fdtd_solver_id == 0) {
- // CFL time step Yee solver
- deltat = cfl * 1./( std::sqrt(AMREX_D_TERM( 1./(dx[0]*dx[0]),
+ // CFL time step Yee solver
+ deltat = cfl * 1./( std::sqrt(AMREX_D_TERM( 1./(dx[0]*dx[0]),
+ 1./(dx[1]*dx[1]),
+ 1./(dx[2]*dx[2]))) * PhysConst::c );
} else {
- // CFL time step CKC solver
+ // CFL time step CKC solver
#if (BL_SPACEDIM == 3)
- const Real delta = std::min(dx[0],std::min(dx[1],dx[2]));
+ const Real delta = std::min(dx[0],std::min(dx[1],dx[2]));
#elif (BL_SPACEDIM == 2)
- const Real delta = std::min(dx[0],dx[1]);
+ const Real delta = std::min(dx[0],dx[1]);
#endif
- deltat = cfl*delta/PhysConst::c;
+ deltat = cfl*delta/PhysConst::c;
}
dt.resize(0);
dt.resize(max_level+1,deltat);