The documentation says the = operator does interpolation. But not clear on what happens when interpolating discontinuous fields. For example:

```
mesh Th = square(1,1);
fespace Vh0(Th,P0);
fespace Vh1(Th,P1);
Vh0 F0 = x;
Vh1 F1 = F0;
Vh1 F1p, q;
plot(F0,fill=1,value=1,dim=3,wait=1);
plot(F1,fill=1,value=1,dim=3,wait=1);
solve proj(F1p,q) = int2d(Th)(F1p*q)-int2d(Th)(F0*q);
plot(F1p,fill=1,value=1,dim=3);
```

One way to go from P0 to P1 field is using a projection which gives F1p.

But for doing interpolation (F1) the value should be computed on points. I’m wondering how the “=” operator computes the values at diagonals? Here we have two different values and it doesn’t average them. Is it picking randomly one of the values of F0 to assign F1?