Remove particles that are initialized in the EB

[RemiLehe]

The current version of WarpX allows particles to be created in the EB. (In that case, these particles only get removed at the end of the first timestep, which may be confusing to some users. For instance, these particles will appear in the diagnostics associated with t=0, and will deposit charge/current which can be seen in the diagnostic of rho/J.

This PR modifies the particle injection so that these particles can marked (with negative IDs) when they are created, and are removed by the next ReDistribute, which should in general occur before the beginning of the first timestep.

[dpgrote]

This will work but will result in extra work setting up the particles that will be then removed.

Looking at what scrapeParticles is doing, the check for whether a particle is in the EB is a small bit of code and could easily be added in the loop adding the particles so any particles in the EB are immediately rejected. Would this be better?

                int ieb, jeb, keb;
                amrex::Real Web[AMREX_SPACEDIM][2];
                ablastr::particles::compute_weights_nodal(ppos.x, ppos.y, ppos.z, plo, dxi, ieb, jeb, keb, Web);
                amrex::Real phi_value = ablastr::particles::interp_field_nodal(ieb, jeb, keb, Web, phi);
                if (phi_value < 0.0) {
                    p.id() = -1;
                    continue;
                }

Though I realize that it wouldn't save much since it would have to happen after the UpdatePosition anyway.

[roelof-groenewald]

Your comment about removing particles after the UpdatePosition call makes me think there is actually a bug in the code presently (besides giving confusing diagnostics on the first step). If I initialize a thermal plasma everywhere in the domain (including over the EB), but we don't remove those immediately, some particles (inside the EB) with velocities directed across the boundary will make it into the valid domain on the first push. Those particles will therefore lead to an undesired preferential velocity (away from the EB) along the edge of the boundary. I guess all this is to say that particles should be removed before the initial push (and also before UpdatePosition, right?).

[dpgrote]

My comment about UpdatePosition only applies to the flux injection, where particles are loaded onto a plane and UpdatePosition is used to advance the particles a fraction of a time step off of the plane. But your comment would still apply, so it would be better to do the check before the position update, which means doing the check as I suggested. This would also mean that any particles that started from a valid location but then crossed into an EB on its initial particle step would be included in the lost particle diagnostic.

[RemiLehe]

Good points ; we could create a small device function isInsideEB to do this.

[RemiLehe]

From the above discussion, we could implement the following changes in a second PR:

• Replacing scrapeParticles by the code that Dave suggested and which should be more economical
• For flux injection: Checking whether particles are in the EB before the UpdatePosition (https://github.com/ECP-WarpX/WarpX/blob/development/Source/Particles/PhysicalParticleContainer.cpp#L1895)

[roelof-groenewald]

As discussed on Slack, we will merge this now so other PRs can move forward and the issues pointed out here will be fixed in subsequent PRs.