From d2cd531cfd1644c15bc9a9852e07931fc5f86810 Mon Sep 17 00:00:00 2001 From: Severin Diederichs <65728274+SeverinDiederichs@users.noreply.github.com> Date: Wed, 17 Mar 2021 19:25:59 +0100 Subject: [PATCH] add option to not deposit Jx and Jy of the beam (#431) --- doc/source/run/parameters.rst | 4 ++++ src/Hipace.H | 2 ++ src/Hipace.cpp | 5 ++++- src/particles/MultiBeam.H | 4 +++- src/particles/MultiBeam.cpp | 6 ++++-- src/particles/deposition/BeamDepositCurrent.H | 4 +++- .../deposition/BeamDepositCurrent.cpp | 19 ++++++++++--------- .../deposition/BeamDepositCurrentInner.H | 18 +++++++++++------- 8 files changed, 41 insertions(+), 21 deletions(-) diff --git a/doc/source/run/parameters.rst b/doc/source/run/parameters.rst index 48e74492ce..30342a9091 100644 --- a/doc/source/run/parameters.rst +++ b/doc/source/run/parameters.rst @@ -47,6 +47,10 @@ General parameters requires | `hipace.beam_injection_cr = 8`. +* ``hipace.do_beam_jx_jy_deposition`` (`bool`) optional (default `1`) + Using the default, the beam deposits all currents `Jx`, `Jy`, `Jz`. Using + `hipace.do_beam_jx_jy_deposition = 0` disables the transverse current deposition of the beams. + Predictor-corrector loop parameters ----------------------------------- diff --git a/src/Hipace.H b/src/Hipace.H index a8ba083010..7fe1ea76bb 100644 --- a/src/Hipace.H +++ b/src/Hipace.H @@ -209,6 +209,8 @@ public: /** Mixing factor between the transverse B field iterations in the predictor corrector loop */ static amrex::Real m_predcorr_B_mixing_factor; + /** Whether the beams deposit Jx and Jy */ + static bool m_do_beam_jx_jy_deposition; /** Whether to call amrex::Gpu::synchronize() around all profiler region */ static bool m_do_device_synchronize; /** How much the box is coarsened for beam injection, to avoid exceeding max int in cell count. diff --git a/src/Hipace.cpp b/src/Hipace.cpp index e8aa9c4762..73074e259c 100644 --- a/src/Hipace.cpp +++ b/src/Hipace.cpp @@ -33,6 +33,7 @@ int Hipace::m_depos_order_z = 0; amrex::Real Hipace::m_predcorr_B_error_tolerance = 4e-2; int Hipace::m_predcorr_max_iterations = 30; amrex::Real Hipace::m_predcorr_B_mixing_factor = 0.05; +bool Hipace::m_do_beam_jx_jy_deposition = true; bool Hipace::m_do_device_synchronize = false; int Hipace::m_beam_injection_cr = 1; amrex::Real Hipace::m_external_ExmBy_slope = 0.; @@ -85,6 +86,7 @@ Hipace::Hipace () : AMREX_ALWAYS_ASSERT_WITH_MESSAGE(m_numprocs_x*m_numprocs_y*m_numprocs_z == amrex::ParallelDescriptor::NProcs(), "Check hipace.numprocs_x and hipace.numprocs_y"); + pph.query("do_beam_jx_jy_deposition", m_do_beam_jx_jy_deposition); pph.query("do_device_synchronize", m_do_device_synchronize); pph.query("external_ExmBy_slope", m_external_ExmBy_slope); pph.query("external_Ez_slope", m_external_Ez_slope); @@ -388,7 +390,8 @@ Hipace::SolveOneSlice (int islice, int lev, const int ibox, m_fields.SolvePoissonExmByAndEypBx(Geom(lev), m_comm_xy, lev); m_grid_current.DepositCurrentSlice(m_fields, geom[lev], lev, islice); - m_multi_beam.DepositCurrentSlice(m_fields, geom[lev], lev, islice, bx, bins, m_box_sorters, ibox); + m_multi_beam.DepositCurrentSlice(m_fields, geom[lev], lev, islice, bx, bins, m_box_sorters, + ibox, m_do_beam_jx_jy_deposition); j_slice.FillBoundary(Geom(lev).periodicity()); diff --git a/src/particles/MultiBeam.H b/src/particles/MultiBeam.H index fe4e94a6a6..0bc9cf7f7d 100644 --- a/src/particles/MultiBeam.H +++ b/src/particles/MultiBeam.H @@ -29,11 +29,13 @@ public: * \param[in] bins Vector (over species) of particles sorted by slices * \param[in] a_box_sorter_vec Vector (over species) of particles sorted by box * \param[in] ibox index of the current box + * \param[in] do_beam_jx_jy_deposition whether the beam deposits Jx and Jy */ void DepositCurrentSlice ( Fields& fields, const amrex::Geometry& geom, const int lev, int islice, const amrex::Box bx, amrex::Vector> bins, - const amrex::Vector& a_box_sorter_vec, const int ibox); + const amrex::Vector& a_box_sorter_vec, const int ibox, + const bool do_beam_jx_jy_deposition); /** Loop over all beam species and build and return the indices of particles sorted per slice * \param[in] lev MR level diff --git a/src/particles/MultiBeam.cpp b/src/particles/MultiBeam.cpp index 6252e26b81..afa45526bb 100644 --- a/src/particles/MultiBeam.cpp +++ b/src/particles/MultiBeam.cpp @@ -27,12 +27,14 @@ void MultiBeam::DepositCurrentSlice ( Fields& fields, const amrex::Geometry& geom, const int lev, int islice, const amrex::Box bx, amrex::Vector> bins, - const amrex::Vector& a_box_sorter_vec, const int ibox) + const amrex::Vector& a_box_sorter_vec, const int ibox, + const bool do_beam_jx_jy_deposition) { for (int i=0; i& bins); + amrex::DenseBins& bins, + const bool do_beam_jx_jy_deposition); #endif // BEAMDEPOSITCURRENT_H_ diff --git a/src/particles/deposition/BeamDepositCurrent.cpp b/src/particles/deposition/BeamDepositCurrent.cpp index fbb5045a0b..7bfbad5ac7 100644 --- a/src/particles/deposition/BeamDepositCurrent.cpp +++ b/src/particles/deposition/BeamDepositCurrent.cpp @@ -11,7 +11,8 @@ void DepositCurrentSlice (BeamParticleContainer& beam, Fields& fields, amrex::Geometry const& gm, int const lev ,const int islice, const amrex::Box bx, int const offset, - amrex::DenseBins& bins) + amrex::DenseBins& bins, + const bool do_beam_jx_jy_deposition) { HIPACE_PROFILE("DepositCurrentSlice_BeamParticleContainer()"); // Extract properties associated with physical size of the box @@ -51,17 +52,17 @@ DepositCurrentSlice (BeamParticleContainer& beam, Fields& fields, amrex::Geometr // Call deposition function in each box if (Hipace::m_depos_order_xy == 0){ - doDepositionShapeN<0, 0>( beam, jx_fab, jy_fab, jz_fab, - dx, xyzmin, lo, q, islice_local, bins, offset); + doDepositionShapeN<0, 0>( beam, jx_fab, jy_fab, jz_fab, dx, xyzmin, lo, q, islice_local, + bins, offset, do_beam_jx_jy_deposition); } else if (Hipace::m_depos_order_xy == 1){ - doDepositionShapeN<1, 0>( beam, jx_fab, jy_fab, jz_fab, - dx, xyzmin, lo, q, islice_local, bins, offset); + doDepositionShapeN<1, 0>( beam, jx_fab, jy_fab, jz_fab, dx, xyzmin, lo, q, islice_local, + bins, offset, do_beam_jx_jy_deposition); } else if (Hipace::m_depos_order_xy == 2){ - doDepositionShapeN<2, 0>( beam, jx_fab, jy_fab, jz_fab, - dx, xyzmin, lo, q, islice_local, bins, offset); + doDepositionShapeN<2, 0>( beam, jx_fab, jy_fab, jz_fab, dx, xyzmin, lo, q, islice_local, + bins, offset, do_beam_jx_jy_deposition); } else if (Hipace::m_depos_order_xy == 3){ - doDepositionShapeN<3, 0>( beam, jx_fab, jy_fab, jz_fab, - dx, xyzmin, lo, q, islice_local, bins, offset); + doDepositionShapeN<3, 0>( beam, jx_fab, jy_fab, jz_fab, dx, xyzmin, lo, q, islice_local, + bins, offset, do_beam_jx_jy_deposition); } else { amrex::Abort("unknown deposition order"); } diff --git a/src/particles/deposition/BeamDepositCurrentInner.H b/src/particles/deposition/BeamDepositCurrentInner.H index 86b377b905..929f178235 100644 --- a/src/particles/deposition/BeamDepositCurrentInner.H +++ b/src/particles/deposition/BeamDepositCurrentInner.H @@ -35,6 +35,7 @@ * \param[in] islice Particles in slice islice will deposit * \param[in] bins Indices of particles arranged per slices. * \param[in] box_offset offset to particles on this box. + * \param[in] do_beam_jx_jy_deposition whether the beams deposit Jx and Jy */ template void doDepositionShapeN (const BeamParticleContainer& ptile, @@ -47,7 +48,8 @@ void doDepositionShapeN (const BeamParticleContainer& ptile, amrex::Real const q, int islice, amrex::DenseBins& bins, - int box_offset) + int box_offset, + const bool do_beam_jx_jy_deposition) { using namespace amrex::literals; @@ -153,12 +155,14 @@ void doDepositionShapeN (const BeamParticleContainer& ptile, for (int iz=0; iz<=depos_order_z; iz++){ for (int iy=0; iy<=depos_order_xy; iy++){ for (int ix=0; ix<=depos_order_xy; ix++){ - amrex::Gpu::Atomic::Add( - &jx_arr(lo.x+j_cell+ix, lo.y+k_cell+iy, z_slice), - sx_cell[ix]*sy_cell[iy]*sz_cell[iz]*wqx); - amrex::Gpu::Atomic::Add( - &jy_arr(lo.x+j_cell+ix, lo.y+k_cell+iy, z_slice), - sx_cell[ix]*sy_cell[iy]*sz_cell[iz]*wqy); + if (do_beam_jx_jy_deposition) { + amrex::Gpu::Atomic::Add( + &jx_arr(lo.x+j_cell+ix, lo.y+k_cell+iy, z_slice), + sx_cell[ix]*sy_cell[iy]*sz_cell[iz]*wqx); + amrex::Gpu::Atomic::Add( + &jy_arr(lo.x+j_cell+ix, lo.y+k_cell+iy, z_slice), + sx_cell[ix]*sy_cell[iy]*sz_cell[iz]*wqy); + } amrex::Gpu::Atomic::Add( &jz_arr(lo.x+j_cell+ix, lo.y+k_cell+iy, z_slice), sx_cell[ix]*sy_cell[iy]*sz_cell[iz]*wqz);