Operations between objects on the boundary and not

[zjwegert]

In #176 I encountered the following issue when trying to do operations between an FEFunction defined via a BoundaryTriangulation and another that is on the background model. Here is an MWE:

using Gridap
using GridapDistributed
using PartitionedArrays

np = (2,2)
ranks = with_debug() do distribute
  distribute(LinearIndices((prod(np),)))
end

domain = (0,1,0,1)
partition = (5,5)
model = CartesianDiscreteModel(ranks,np,domain,partition)
Γ = BoundaryTriangulation(model)
V1 = FESpace(Γ,ReferenceFE(lagrangian,Float64,1))
V2 = FESpace(model,ReferenceFE(lagrangian,Float64,1))
uh1 = zero(V1);
uh2 = zero(V2);
uh1*uh2

This produces an unreachable error because _to_common_domain fails to find a common domain even though it should. This appears to be because of the following: for DistributedTriangulations resulting from BoundaryTriangulation, we're in Case 4 of add_ghost_cells, so we create a new create a new triangulation with the ghost cells. The issue then appears to be that the glue information is missing, so both cases of

is_change_possible(_trian_a,_trian_b)
is_change_possible(_trian_b,_trian_a)

return false.

[JordiManyer]

Your assessment is correct, I am thinking on how to solve it in the most elegant way possible. This ghosted/non-ghosted issue has been popping up quite a bit, so I want to actually think of a solution.

[JordiManyer]

Actually this is repeated from #170

[oriolcg]

We started having a look at this in `mixdims' branch a while back with @fverdugo, but never managed to finish. Maybe there are things there that could be useful to have as starting point.

[JordiManyer]

Its not hard to solve per se, its just hard to solve elegantly I think.

We are currently returning the correct BodyFittedTriangulation, it's just a matter of communicating the glue (which can be done easily).

What I am pondering is whether this is what we want. For Boundary/Skeleton trians, we are actually keeping a view of the unghosted triangulation. So in theory we just need to take the parent and we are fine. The issue is that, right now, there is no way of knowing whether the you are in this case or not without communicating.
So what I've been thinking is whether we should add some kind of trait/metadata that carries this info in the cases where this is true. Or maybe for the sake of generality we just do the communication and thats it. In the end the price is payed once (in most cases), so it might make sense to do something generic.

[fverdugo]

As far as I remember, the toughest problem was in the operations on cell fields for distributed mix dims interpolation spaces.

[JordiManyer]

OK so quick poll: When you call

trian = Boundary(WithGhost(),model)

what would you expect to get?

• The local part of the global physical boundary, or
• The local physical boundary, i.e with the interface between processors included?

The difference between the two being the presence of the processor interface faces or not.
Right now, in the master branch we are getting the latter, but I would almost call this a bug?