Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

Minor refactoring of refined injection for AddPlasmaFlux #5274

Merged
merged 5 commits into from
Sep 17, 2024
Merged

Minor refactoring of refined injection for AddPlasmaFlux #5274

Changes from all commits
Commits
Show all changes
5 commits
Select commit Hold shift + click to select a range
89c684f
Minor refactoring of refined injection for AddPlasmaFlux
RemiLehe Sep 15, 2024
3c00605
Additional simplification
RemiLehe Sep 16, 2024
f16b01c
Minor fixes
RemiLehe Sep 16, 2024
bd086a1
Merge branch 'development' into simplify_fine_injectionwq
RemiLehe Sep 16, 2024
c6d8ac9
Remove unneeded variable lrrfac
RemiLehe Sep 17, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
37 changes: 12 additions & 25 deletions Source/Particles/PhysicalParticleContainer.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -1030,7 +1030,6 @@ PhysicalParticleContainer::AddPlasma (PlasmaInjector const& plasma_injector, int
Gpu::DeviceVector<amrex::Long> counts(overlap_box.numPts(), 0);
Gpu::DeviceVector<amrex::Long> offset(overlap_box.numPts());
auto *pcounts = counts.data();
const amrex::IntVect lrrfac = rrfac;
Box fine_overlap_box; // default Box is NOT ok().
if (refine_injection) {
fine_overlap_box = overlap_box & amrex::shift(fine_injection_box, -shifted);
Expand All @@ -1048,7 +1047,7 @@ PhysicalParticleContainer::AddPlasma (PlasmaInjector const& plasma_injector, int
{
auto index = overlap_box.index(iv);
const amrex::Long r = (fine_overlap_box.ok() && fine_overlap_box.contains(iv))?
(AMREX_D_TERM(lrrfac[0],*lrrfac[1],*lrrfac[2])) : (1);
(AMREX_D_TERM(rrfac[0],*rrfac[1],*rrfac[2])) : (1);
pcounts[index] = num_ppc*r;
// update pcount by checking if cell-corners or cell-center
// has non-zero density
Expand Down Expand Up @@ -1154,8 +1153,8 @@ PhysicalParticleContainer::AddPlasma (PlasmaInjector const& plasma_injector, int
long ip = poffset[index] + i_part;
pa_idcpu[ip] = amrex::SetParticleIDandCPU(pid+ip, cpuid);
const XDim3 r = (fine_overlap_box.ok() && fine_overlap_box.contains(iv)) ?
// In the refined injection region: use refinement ratio `lrrfac`
inj_pos->getPositionUnitBox(i_part, lrrfac, engine) :
// In the refined injection region: use refinement ratio `rrfac`
inj_pos->getPositionUnitBox(i_part, rrfac, engine) :
// Otherwise: use 1 as the refinement ratio
inj_pos->getPositionUnitBox(i_part, amrex::IntVect::TheUnitVector(), engine);
auto pos = getCellCoords(overlap_corner, dx, r, iv);
Expand Down Expand Up @@ -1441,7 +1440,6 @@ PhysicalParticleContainer::AddPlasmaFlux (PlasmaInjector const& plasma_injector,
Gpu::DeviceVector<int> counts(overlap_box.numPts(), 0);
Gpu::DeviceVector<int> offset(overlap_box.numPts());
auto *pcounts = counts.data();
const amrex::IntVect lrrfac = rrfac;
const int flux_normal_axis = plasma_injector.flux_normal_axis;
Box fine_overlap_box; // default Box is NOT ok().
if (refine_injection) {
Expand All @@ -1450,22 +1448,21 @@ PhysicalParticleContainer::AddPlasmaFlux (PlasmaInjector const& plasma_injector,
amrex::ParallelForRNG(overlap_box, [=] AMREX_GPU_DEVICE (int i, int j, int k, amrex::RandomEngine const& engine) noexcept
{
const IntVect iv(AMREX_D_DECL(i, j, k));
amrex::ignore_unused(j,k);

auto lo = getCellCoords(overlap_corner, dx, {0._rt, 0._rt, 0._rt}, iv);
auto hi = getCellCoords(overlap_corner, dx, {1._rt, 1._rt, 1._rt}, iv);

if (flux_pos->overlapsWith(lo, hi))
{
auto index = overlap_box.index(iv);
int r;
int r = 1;
if (fine_overlap_box.ok() && fine_overlap_box.contains(iv)) {
r = compute_area_weights(lrrfac, flux_normal_axis);
} else {
r = 1;
r = compute_area_weights(rrfac, flux_normal_axis);
}
const int num_ppc_int = static_cast<int>(num_ppc_real*r + amrex::Random(engine));
pcounts[index] = num_ppc_int;
}
amrex::ignore_unused(j,k);
});

// Max number of new particles. All of them are created,
Expand Down Expand Up @@ -1532,24 +1529,14 @@ PhysicalParticleContainer::AddPlasmaFlux (PlasmaInjector const& plasma_injector,
[=] AMREX_GPU_DEVICE (int i, int j, int k, amrex::RandomEngine const& engine) noexcept
{
const IntVect iv = IntVect(AMREX_D_DECL(i, j, k));
amrex::ignore_unused(j,k);
const auto index = overlap_box.index(iv);

Real scale_fac = compute_scale_fac_area(dx, num_ppc_real, flux_normal_axis);

auto lo = getCellCoords(overlap_corner, dx, {0._rt, 0._rt, 0._rt}, iv);
auto hi = getCellCoords(overlap_corner, dx, {1._rt, 1._rt, 1._rt}, iv);

if (flux_pos->overlapsWith(lo, hi))
Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This if condition is not key here, since the pcounts[index] is 0 anyway for cells that do not overlap.

{
int r;
if (fine_overlap_box.ok() && fine_overlap_box.contains(iv)) {
r = compute_area_weights(lrrfac, flux_normal_axis);
} else {
r = 1;
}
scale_fac /= r;
if (fine_overlap_box.ok() && fine_overlap_box.contains(iv)) {
scale_fac /= compute_area_weights(rrfac, flux_normal_axis);
}
amrex::ignore_unused(j,k);

for (int i_part = 0; i_part < pcounts[index]; ++i_part)
{
Expand All @@ -1558,8 +1545,8 @@ PhysicalParticleContainer::AddPlasmaFlux (PlasmaInjector const& plasma_injector,

// This assumes the flux_pos is of type InjectorPositionRandomPlane
const XDim3 r = (fine_overlap_box.ok() && fine_overlap_box.contains(iv)) ?
// In the refined injection region: use refinement ratio `lrrfac`
flux_pos->getPositionUnitBox(i_part, lrrfac, engine) :
// In the refined injection region: use refinement ratio `rrfac`
flux_pos->getPositionUnitBox(i_part, rrfac, engine) :
// Otherwise: use 1 as the refinement ratio
flux_pos->getPositionUnitBox(i_part, amrex::IntVect::TheUnitVector(), engine);
auto pos = getCellCoords(overlap_corner, dx, r, iv);
Expand Down
Loading