Hello. I am trying to implement a 1D Poisson problem with weak form a(u,v)=L(v), where the right hand side is given by L(v)=-v'(0). From searching for similar questions in this forum, I arrived at the following idea to implement the right hand side, where the idea is to approximate the dirac delta.
h = 0.8E-3
delta = Expression('h*pow(pi,-1)*pow((pow(x[0],2) + pow(h,2)),-1)',domain=mesh,degree = 4, h = h)
L = -delta*grad(v)[0]*dx
However, I would have to somehow select a good value for h in this case, which is a bit unsatisfactory. Is there a better/prefered/easier way to implement such a right hand side in Fenics?
Before implementing, it’s worth stepping back and noting a theoretical problem. This is the Poisson problem with a Dirac delta derivative \delta' at the source term, but \delta'\notin H^{-1}, so you end up outside the usual variational framework of considering test functions to be in H^1, because \langle\delta',v\rangle is not well-defined for all v\in H^1. You can also think about it in the discrete setting: If x=0 were a node of a continuous finite element space, v'(0) would not be well-defined. Your exact solution will have a discontinuity at x=0, so it won’t be in H^1.
You might try using an analytical singularity-removal technique, analogous to this, but simpler in 1D. You can solve exactly for the solution with just the Dirac delta derivative as a source term, to get something analogous to a Green’s function, but for the Dirac delta derivative:
where H(\cdot) is the Heaviside function and C_1 and C_2 are integration constants that can be chosen to satisfy homogeneous boundary conditions at the ends of the interval. Then you can reformulate your problem with
u~\to~u_\text{reg} + G\text{ ,}
where the unknown becomes the regular part of the solution, u_\text{reg}, subject to inhomogenous boundary conditions:
where f_\text{reg} is the regular part of the source term (and, if understood in a generalized sense, can be though of as a functional including any inhomogeneous boundary conditions).