On Wed, 2 Mar 2005, John Peterson wrote:
> Currently, the Kelly Error Estimator ignores contributions to the
> residual coming from element sides on the boundary which have Neumann
> (flux) type boundary conditions such as
> grad(u).n = g on \Gamma_N
> I think it is reasonable to assume this error will be small, since the
> flux condition is enforced weakly during the formation of the element
> stiffness matrices. In the case of natural (flux=0) boundary conditions,
> it is common to not even add these terms to the stiffness matrix at all
> since they are assumed to be zero.
You know, although it seems natural to assume the flux error will be
small, I'm not sure how rapid the convergence rate is - in fact I
wouldn't be surprised if there's some Trace Theorem nastiness that
loses us half a power of h. In any case, I don't think there's any
reason to assume the boundary flux error will be any better than the
interior flux jumps.
> My question is: would it be worth computing the boundary edge integrals
> \int |g - grad(u_h).n|^2 dS
> on the flux boundaries? My hunch is that in some applications, the
> estimator should actually be reporting a lot more error for the
> elements on the boundary than it currently is. In fact, I have seen
> actual evidence of this in a 2D tumor-growth model with an exponential
> layer near the boundary. The row of elements closest to the boundary
> are not flagged for refinement even though several adjacent layers of
> elements are. What are your thoughts?
That sounds like a problem to me.
> The main implementation difficulty would be coming up with a general
> way to get the BCs and BC information into the error estimator without
> duplicating existing code in the matrix assembly routine. There are
> probably other issues as well. Any thoughts would be appreciated.
This is probably something that could be addressed when we're
designing the new system assembly interface. Currently libmesh
doesn't really know about boundary conditions - only the user's
assembly function knows. Ideally the assembly function would be
handed boundary conditions from the system along with the rest of the
problem description, and the user would assign those boundary
conditions to the system in main(). That could make it easier to
enforce Dirichlet conditions exactly instead of with penalty methods,