handle in-patch invariants when converting between different precisions

[germasch]

So I'm pretty sure I hit cases where a particle was validly inside a given patch in double precision, but when those particle got converted to single precision for some operation (like sorting), that invariant was violated and then leads to an assertion failure (or worse).

[wrfox]

Kai, In the Fortran version, PIC_pe* is called immediately before sort, to make sure all particles are on the correct proc. I am sure this is not optimal, but is guaranteed to be correct. If switching from double to single can put the simulation into an inconsistent state, I guess that the C version analogy of PIC_pe* has to occur every time the precision is downgraded? Will

…

On Feb 6, 2018, at 7:41 AM, Kai Germaschewski ***@***.***> wrote: So I'm pretty sure I hit cases where a particle was validly inside a given patch in double precision, but when those particle got converted to single precision for some operation (like sorting), that invariant was violated and then leads to an assertion failure (or worse). — You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHub, or mute the thread.

[germasch]

Yeah, so I think it helps to have proper invariants, and I've moved into that direction: Ideally, the invariant would just be: All particles in the particle array are within bounds (which needs a definition, too).

In particular, whenever a particle gets added (during injection, initialization, ?), it needs to be within bounds. However, the particle push breaks the invariant temporarily, until the boundary exchange is performed. This should be kinda acceptable, though not ideal, since it's a bit more complex than I'd like to have it. One option might be to move particles that leave the local patch out out of the particle list right in the pusher, so then the invariant would be satisfied at all times.

The particular issues at hand, though, I think isn't as tough. Patch ownership is currently defined as, say the domain goes from [0:40[, divided into 4 subdomains / patches in that direction:
[0:10[, [10:20[, [20:30[, [30:40[. Patch positions are actually stored patch-local, so local positions would be [0:10[ for every patch. However, if one does not have periodic b.c., the last patch should actually become [30:40], or [0:10] locally. That's kind of a pain, since the cell index isn't then just simply floor(x / dx), a particle exactly at the right limit would be considered in a non-existing cell. This is still not handled well.

The particular problem here, though, is not that bad, it doesn't require actual communication to solve, I think. A particle at 9.999999999 in double precision turns into position 10.0 in single precision, which would be considered out of bounds. But I think just nudging it back to 9.99999 would work well enough.

On the other hand, I think the best solution might be to just have the invariant condition to allow the closed interval [0:10] everywhere, but assign a particle at the right limit to be still at the cell index m-1, that is, in the cell [9:10[. Interpolation etc will work just fine that way, as long as that cell index is always calculated consistently that way.

I don't completely like that because of the ambiguity (a particle at 20.0 could legally part of patch [10:20] or [20:30]), but since some kind of special case, at least at the top boundary w/o periodic b.c. is unavoidable, that may be the cleanest way.

[wrfox]

Hi Kai, Okay, so it sounds like particle’s cell (and host patch) is changing if the position is truncated to single first. Effectively they go into a ghost cell on a patch. I have to say, though, what you write is happening should not make a difference on the physics side for the particle pusher The interpolation I recall is constructed s.t. it barely matters if a particle exactly lands on one side of a cell boundary or another - because the relevant interpolation coefficients are going to zero as a particle approaches the boundary. Therefore, on the physics side, at least, it seems to me that a particle on the knife edge of a cell should be allowed to be in either cell. So this raises the question to me of why is this check required? (However, in contrast - in the collision routine, we cannot have particles in ghost cells going in, because then the pairing will break - the particles need to be in their true host proc. This is why in the fortran version, we call a PIC_pe* first to make sure no particles are in ghost cells. Will

…

On Feb 6, 2018, at 1:40 PM, Kai Germaschewski ***@***.***> wrote: Yeah, so I think it helps to have proper invariants, and I've moved into that direction: Ideally, the invariant would just be: All particles in the particle array are within bounds (which needs a definition, too). In particular, whenever a particle gets added (during injection, initialization, ?), it needs to be within bounds. However, the particle push breaks the invariant temporarily, until the boundary exchange is performed. This should be kinda acceptable, though not ideal, since it's a bit more complex than I'd like to have it. One option might be to move particles that leave the local patch out out of the particle list right in the pusher, so then the invariant would be satisfied at all times. The particular issues at hand, though, I think isn't as tough. Patch ownership is currently defined as, say the domain goes from [0:40[, divided into 4 subdomains / patches in that direction: [0:10[, [10:20[, [20:30[, [30:40[. Patch positions are actually stored patch-local, so local positions would be [0:10[ for every patch. However, if one does not have periodic b.c., the last patch should actually become [30:40], or [0:10] locally. That's kind of a pain, since the cell index isn't then just simply floor(x / dx), a particle exactly at the right limit would be considered in a non-existing cell. This is still not handled well. The particular problem here, though, is not that bad, it doesn't require actual communication to solve, I think. A particle at 9.999999999 in double precision turns into position 10.0 in single precision, which would be considered out of bounds. But I think just nudging it back to 9.99999 would work well enough. On the other hand, I think the best solution might be to just have the invariant condition to allow the closed interval [0:10] everywhere, but assign a particle at the right limit to be still at the cell index m-1, that is, in the cell [9:10[. Interpolation etc will work just fine that way, as long as that cell index is always calculated consistently that way. I don't completely like that because of the ambiguity (a particle at 20.0 could legally part of patch [10:20] or [20:30]), but since some kind of special case, at least at the top boundary w/o periodic b.c. is unavoidable, that may be the cleanest way. — You are receiving this because you commented. Reply to this email directly, view it on GitHub, or mute the thread.