Streamlines using the Elmer method

[fb77 (François Bouchut)]

A method is proposed in

Uzawa scheme for stokes bingham flow

In that paper the flow is taken incompressible, whereas in your code you use compressible. This is why the solution is different (but looks correct). You can modify easily your code to do the incompressible case, by using the extra variables p,q as solve bing([ux,uy,p],[vx,vy,q]) and adding -int2d(Th)(p*Divergence(vx, vy)+q*Divergence(ux, uy)) But then it happens that the scheme converges extremely slowly (you need thousands of iterations). In order to improve it you should use the Augmente…

from “Another point is that your computation of the stream function…”

The important thing is the \varepsilon.
If your u is P1, you should take the stream function in P2.
If your u is P2, you should take the stream function in P3.