Problems with error calculation for solution with singularity

NGSolve
1

Hi all,

I’m writing my bachelor’s thesis on AFEM at the moment; my focus is on comparing the performance of a residual based error estimation approach and the ZZ estimation approach for a simple Poisson problem. I’m currently trying to implement both into NGSolve, I thought I had a working version testing it on a simple problem with a known solution (poisson problem with u = \sin(\pi x)\sin(\pi y) on an L-shaped domain.

This doesn’t yield super interesting results, so I’m trying to implement it with a solution that should have a singularity. My supervisor suggested I try u = r^{2/3} \sin(\frac{2}{3}\theta)(1-x^2)(1-y^2) (the Cartesian terms are there to meet Dirichlet boundary conditions).

I’m trying to implement this problem into my solver, but I get unexpected behaviour; the total error does not decrease per AFEM iteration, even though the error estimation and mesh refinement steps seem to be running fine. (This weird behaviour of course does suggest an error, I just can’t seem to find it…).

I’ve attached my code (tried to summarise a bit, but not quite sure where my error is…), any tips would be greatly appreciated. My hunch is that it’s got something to do with my integration steps to get the exact error over each element as this could be inaccurate near a singularity?

Solver_Share.ipynb (34.3 KB)

2

Hey! One possible reason your total error isn’t decreasing could be that the singular solution:

u=r^{2/3}\sin(\frac{2}{3}\theta)(1-x^2)(1-y^2)

doesn’t actually induce a true singularity on the L-shaped domain — the (1-x^2)(1-y^2) term smooths it out too much. So AFEM isn’t refining near the re-entrant corner like it should. I’d suggest switching back to a smooth solution like \sin(\pi x)\sin(\pi y) if you want the estimators to behave better. Also, you probably don’t need to worry about integration accuracy for the exact error since the FEM solution should already “know” how to handle singularities inherently.

3

Hmmm, ok, I’ll have a look into that. I have managed to get the error to decrease with AFEM iterations now (by switching to atan2(y,x) from tan(y/x) to get \theta). However the convergence is still very slow…

Any ideas why that might be?

4

This might be not an issue at all, but you may want to double-check what the branch of \operatorname{arctan} the function atan2 uses.

1 Like
5

This was it; I needed to have \theta \in[0,\frac{3\pi}{2}]. Thank you for the help!

1 Like
6

Glad to help! I had this same exact issue years ago when solving a homework!