Point source not working

Hi,

I want to use a Dirac function in my code and I want to just print the result before using it. I use an 1D mesh and I want to apply it on the node x = 0 but I am getting an error after the projection. Here is my code


mesh = IntervalMesh ( 100 , 0., 1.0 )

class Delta(UserExpression):
    def __init__(self, eps, x0, **kwargs):
        self.eps = eps
        self.x0 = x0
        UserExpression.__init__(self, **kwargs)
    def eval(self, value, x):
        eps = self.eps
        value = eps/pi/(np.linalg.norm(x-self.x0)**2 + eps**2)

    def value_shape(self): return (1, )


delta = Delta(eps=1E-4, x0=0., degree=1)

V = FunctionSpace(mesh, "CG", 1)
u  = project(delta, V)
File("u.pvd") << u

Remove this. Your function is scalar, not a vector of size 1.

Or project into a vector function space of shape (1,).

Although it works, in paraview I can’t see a dirac delta distribution

You’re overwriting the value variable. You should be populating it. In essence:

value[0] = eps/pi/(np.linalg.norm(x-self.x0)**2 + eps**2)

which upon setting x0 = 0.5 gives me

And, if you’re not already aware, the PointSource class offers the integration of Dirac delta function multiplied by a test function.

1 Like

Thank you very much!