Hello everybody,

I have implemented a Dirac Delta UserExpression on dolfin legacy, however, it proves to be quite a slow implementation. I know that the points that I’ll be applying my Dirac Delta to are located over the boundary of my mesh. Is there a way to define the class Delta so that it only works in the boundary, making my code a little faster?

My code:

```
mesh = UnitSquareMesh(100,100)
V = VectorFunctionSpace(mesh, 'Lagrange', degree=1)
class Delta(UserExpression):
def __init__(self, eps, x0, scale_x, scale_y, **kwargs):
self.eps = eps
self.x0 = x0
self.scale_x = scale_x
self.scale_y = scale_y
UserExpression.__init__(self, **kwargs)
def eval(self, values, x):
eps = self.eps
dist_squaredx = np.linalg.norm(x[0] - self.x0[0])**2 + eps**2
dist_squaredy = np.linalg.norm(x[1] - self.x0[1])**2 + eps**2
values[0] = self.scale_x * eps * abs(x[0] - self.x0[0]) / (pi * dist_squaredx)
values[1] = self.scale_y * eps * abs(x[1] - self.x0[1]) / (pi * dist_squaredy)
def value_shape(self):
return (2, )
points = [[0,0] , [1,0] , [0.5,0]]
forces = [[0,-1] , [0,-1], [0,-1]]
# Create a list of Delta expressions for each point and force
deltas = [Delta(eps=1E-4, x0=point, scale_x=force[0]*1000, scale_y=force[1]*1000, degree=5)
for point, force in zip(points, forces)]
```

Thanks in advance