Comments on number of ghost cells

[beck-llr]

This discussion was initiated on the chat by psalz:
I have a bit of a technical question regarding the implementation of the stencil: It is my understanding that each patch has some ghost cells in each dimension that contain data from the neighboring patches. On the "low" side in each dimension there's 2 ghost cells, and on the "high" side there's 2 + 1 + isDual(D). However, when looking at the code for exchanging those values (https://github.com/SmileiPIC/Smilei/blob/d8e1d8595693c1633ab737bb8e0c5364d8f0a1a8/src/Patch/SyncVectorPatch.cpp#L1551-L1564), only oversize (=2) cells are exchanged in each dimension. For the "low" side that is fine, but on the "high" side, only the last two cells are exchanged, i.e., in a field of size N, cells N-1 and N-2. This means that cells N-3 (and if isDual(D), N-4) are NOT exchanged. I looked at the values in a couple of patches and it seems like this works out fine between inner boundaries (N-3 and N-4 still contain the correct values of the neighboring patch), but it doesn't for wrapping around to the other side of the domain when using periodic boundary conditions. In other words, when using periodic boundary conditions, cells 2 and 3 in the first patch do NOT contain the same values as cells N-4 and N-3 of the last patch on the other side of the domain. I tried changing this to do the full oversize + 1 + isDual(D) exchange, which then results in a point-wise error in field values in the order of ~10^-6 and an error in the order of ~10^-4 for Utot, Uelm and Ukin. Could you maybe shed some light on what is going on here?

[beck-llr]

I have several remarks.

When discussing boundaries, please refer to Xmin, Xmax, Ymin, Ymax and so on for clarity because high and low i s ambiguous.

There are the same number of ghost cells in each direction and for all patches. That number is the oversize. There is probably something that you have misunderstand. Within these cells, primal fields have oversize grid points and dual fields have oversize+1 grid points.

When dealing with field exchange you can not just give a cell number because it is ambiguous as well. You need to consider grid points because, as you noticed, quantities are either primal or dual and are not treated equally.

You also have to remember that ghost cells overlap the real cells of the neighbour and not each other. Meaning that with oversize=2 and periodic boundary conditions and for primal grid points going from i=0 to i= N-1 (among which 0,1, N-1 and N-2 are ghosts) you should have F(0) = F(N-5) , F(1)=F(N-4), F(3)=F(N-2) and F(4)=F(N-1). Do we agree on that ?

[psalz]

Thanks for the quick response!

There is probably something that you have misunderstand.

That is very likely, yes.

You also have to remember that ghost cells overlap the real cells of the neighbour and not each other. Meaning that with oversize=2 and periodic boundary conditions and for primal grid points going from i=0 to i= N-1 (among which 0,1, N-1 and N-2 are ghosts) you should have F(0) = F(N-5) , F(1)=F(N-4), F(3)=F(N-2) and F(4)=F(N-1). Do we agree on that ?

I think so, although I think my confusion stems from F(2) and F(N-3) (or on a dual grid F(2), F(3), F(N-3) and F(N-4)) - what are those?

To ease the discussion please have a look at this illustration I've thrown together.

It shows field Bz (dual in X and Y) for a 2D simulation before and after the first call to exchangeB. The grid size is 32x8 cells and there's 4x1 patches. The printed values have been multiplied by 10^6 to make them more readable. I've drawn arrows where transfers are being made, the blue and green boxes indicate matching values.

If we look at the boundary between patch 0,0 and patch 0,1, it seems strange to me that the right (Xmax) ghost cells in patch 0,0 (i.e., F(N-2) and F(N-1)) are NOT the leftmost (Xmin) non-ghost cells (i.e., F(2) and F(3)) of patch 0,1. It seems to me that for all intents and purposes they are ghost cells, as they do in fact contain the same values (indicated by the dark green boxes), even though they are not being exchanged.

They also kind of have to be, because if I have a domain of 32x8 and split it into 4x1 patches (as shown in the picture) I end up with 4 Bz fields of size 14x14 each, or a total of 56x14 cells. If I use a single large patch on the other hand, it is only 38x14 cells. Unless in the 4x1 case there's 4*6 = 24 duplicated cells, for a total of 56 - 24 = 32 normal cells, it would mean that having more patches would result in a larger simulation grid - that doesn't seem right.

So far everything made kind of sense to me, until I got to the part marked in red in the image. After a call to exchangeB, those kind-of-but-not-really ghost cells contain very different values between patches 0,0 and 3,0, when using periodic boundary conditions in X.

What further compounds my confusion is that for summing up J after current projection it seems like all 2*oversize+1+isDual cells are in fact considered ghost cells and are actually being exchanged.

So in conclusion, what are those cells that mostly contain the same data as their neighbor, except when wrapping around the domain for periodic boundary conditions?

[beck-llr]

Very nice illustration. This would indeed be much easier if I could draw on a blackboard.

First I should start by saying that communications are usually minimized in the sense that only necessary information is exchanged. Anything that can be retrieved by mere computation is not exchanged. That explains why less points are exchanged for fields but one more has to be exchanged for currents. Have a look at which grid points can be influenced by a particle of a another patch (currents) and which grid points can not be updated by the local patch solver (fields). These are the points that need information from outside.

it seems strange to me that the right (Xmax) ghost cells in patch 0,0 (i.e., F(N-2) and F(N-1)) are NOT the leftmost (Xmin) non-ghost cells (i.e., F(2) and F(3)) of patch 0,1.
This is dual in X so it is different. Since dual fields extend further away you have patch 0 F(N-1) = patch 1 F(5). Another way to have a look at it is to understand that the last primal grid point (not ghost) i=N-3 is physically the same as the first of the next patch i=3. From there everything overlaps, dual grid is shifted by half a cell towards xmin but has one more grid point so it extends further both in the xmin AND xmax direction.

So in conclusion, what are those cells that mostly contain the same data as their neighbor, except when wrapping around the domain for periodic boundary conditions?
Now I understand your question better and don't have an answer unfortunately. Maybe someone who implemented those boundaries can answer. Normally these points should be correctly evaluated by the solver and do not require communications (as the other dark green areas in your image). Could you please explain how you got these numbers ?

[psalz]

First I should start by saying that communications are usually minimized in the sense that only necessary information is exchanged. Anything that can be retrieved by mere computation is not exchanged. That explains why less points are exchanged for fields but one more has to be exchanged for currents.

Okay thank you, that about matches with what I expected.

Now I understand your question better and don't have an answer unfortunately. Maybe someone who implemented those boundaries can answer. Normally these points should be correctly evaluated by the solver and do not require communications (as the other dark green areas in your image). Could you please explain how you got these numbers ?

The namelist I've used to generate the numbers is included below the image in the gist. The values are printed in VectorPatch::finalizeSyncAndBCFields, right after the call to SyncVectorPatch::finalizeExchangeB. As I've said they're multiplied by 10^6 to make them more legible. I'm still on Smilei 4.8, but it doesn't look like much has changed regarding the boundary exchange.

[psalz]

Anyone else who could chime in on this?