Hello all, I’m new here but not so new to FreeFem++.

I want to implement the convect operator on a periodic function (with a periodic velocity field) in 2d.

mesh Th=square(Nx,Ny,[Lx*x,Ly*y]);

fespace Den(Th,P2,periodic=[[2,y],[4,y],[1,x],[3,x]]);

Den c0=whatever0, c;

fespace Vel(Th,[P1b,P1b],periodic=[[2,y],[4,y],[1,x],[3,x]]);

Vel [u1,u2]=[whatever1,whatever2];

c = convect([u1,u2], -dt, c0);

This does not covect properly through the periodic boundaries (the density on the edges remains static). I am sure that Eulerian advection works in principle but I prefer the Semi-Langrangian route.

The workaround I am considering is:

- Define C0=c0, [U1,U2]=[u1,u2] on a broader mesh that extends beyond the boundaries in all directions. C0 would be periodic over Lx,Ly, same for [U1,U2].
- C = convect([U1,U2], -dt, C0);
- Restrict C to c (ie the original mesh).

Is convect supposed to work on periodic boundaries? If not, is there a more elegant implementation?

Thanks in advance