Applying the Dirac Delta Function w/o using the Fenics PointSource function

Question

Good evening,

I am working on a problem where Dirac has to be applied. I am already aware of the option, assembling the system and then applying the Dirac with the PointSource function on the right hand side vector. My question is, is there another way of applying the Dirac Delta function before or without the system being assembled (before any linear algebra ops take place), eg if I'm using the solve( a == L, u, bcs) function.

Thanks in advance

Comments

Hi, could you pls. tell me what problem are you working on? Thanks.

---

Hi, it's a Poisson problem uxx+uyy=f, in a rectangular domain, Dirichlet BC everywhere but the right side, where Neumann apply, with a point source in the interior of the domain.

Answer 1

Hi, you could use some of the limit definitions of the delta function in the Expression. For example

from dolfin import *

class Delta(Expression):
    def __init__(self, eps):
        self.eps = eps
    def eval(self, values, x):
        eps = self.eps
        values[0] = eps/pi/(x[0]**2 + eps**2) 

mesh = IntervalMesh(1000, -1, 1)
V = FunctionSpace(mesh, 'CG', 1)
bc = DirichletBC(V, Constant(0), 'on_boundary')
u = TrialFunction(V)
v = TestFunction(V)

a = inner(grad(u), grad(v))*dx 
delta = Delta(1E-4)
L = inner(delta, v)*dx

u = Function(V)
if True:
    solve(a == L, u, bc)
else:
    A, b = assemble_system(a, L, bc)
    b.zero()
    PointSource(V, Point(0), 100).apply(b)
    solve(A, u.vector(), b)

plot(delta, mesh=mesh)
plot(u, interactive=True)

Note that in the example I did not tune the parameters of delta or the point source to make the solutions match.

Comments

Thanks for your reply MiroK!