eliminate_internal=True vs False (slightly different results)

I have a question about a different result I’m getting when using eliminate_internal=True vs eliminate_internal=False. As a simple example, say we have Stokes equations but with the mass equation replaced by \varepsilon^{-1} p + \nabla \cdot u = 0. If using VectorL2 with order=2 for velocity and L2 with order=1 for pressure (in an HDG discretization), I expect that if I compute the L^2-norm of \varepsilon^{-1} p_h + \nabla \cdot u_h, once I’ve solved for velocity and pressure, that this norm is (close to) zero. In computing this L^2-norm, with \varepsilon=10^5, this norm is about 8.9e-7 when using eliminate_internal=True and 2.2e-10 when eliminate_internal=False. Is there a simple explanation for this difference? (Maybe extra round-off errors from static condensation?)

Hi Sander,

difficult to say. It certainly depends on the HDG version, I think in your version you use Lagrange parameters for the pressure ?

Do you use iterative solvers - does it depend on the precision ?

What about doing a post-iteration, computing residuals and solve for another update.

How does the quantity depend on epsilon, if you change eps by a factor of 10 larger or smaller. If you see saturation, its likely that you run into floating point precision.

All the best,

Hi Joachim,
Thanks for the response. Yes, I’m using Lagrange multipliers for the pressure. I’m using umfpack as solver. When epsilon = 1, then with eleminate_internal=False I get 5.8e-16 and with eliminate_internal=True I get 1.9e-15. I have attached the code.
stokes_hdg.py (1.9 KB)