Residual-type a posteriori error estimator

Hello,

Assuming that I have the 2D Poisson problem with pure Neumann b.c., is there a way to obtain a residual type estimator for H1
norm, in order to use it in an adaptive algorithm? Namely, I would like to compute the following quantity for every element K
of the mesh,

\[\eta_{1, K}=\left\{h_K^2\left\|f+\Delta u_h\right\|_K^2+\frac{1}{2} \sum_{S \in \partial K \cap S_{int}} h_S\left\|\left[ \nabla u_h\cdot \vec{n}\right]\right\|_S^2+\sum_{S\in \partial K\cap \partial \Omega} h_S\left\|g- \nabla u_h\cdot \vec{n}\right\|_S^2\right\}^{\frac{1}{2}}.\]

Here $f$ and $g$ are the right hand side on the domain and the boundary, respectively, while $S$ is an edge and $S_{int}$ is the set of all  interior edges.

Best regards,
Vex

Hi Vex,
a residual error estimator for a Dirichlet Poisson problem is given by

h = specialcf.mesh_size
n = specialcf.normal(2)

resT = h*h*(Trace(gfu.Operator("hesse"))+f)**2*dx
resE = h*( (grad(gfu)-grad(gfu).Other())*n )**2 *dx(element_vb=BND)
eta1 = Integrate(resT+resE, mesh, element_wise=True)


To adapt to a Neumann problem I think you have to use some kind of indicator function to identify the boundary edges

indBnd = GridFunction(FacetFESpace(mesh,order=0, dirichlet=".*", autoupdate=True), autoupdate=True)
indBnd.Set(1, BND)
resEbnd = indBnd* h*(grad(gfu)*n-g )**2 *dx(element_vb=BND)


indBnd is 1 at the boundary edges and 0 at the inner edges.

Best,
Michael