I have a rather complicated PDE that I'm trying to solve, which requires the evaluation of test and trial functions over cells, internal facets, and points. The first two work fine, via the dx and dS measures. Using the python interface, anything involving a dP returns zero. As an example:
mesh = UnitSquareMesh(2,2)
Q = FunctionSpace(mesh,'CG',1)
phi = TestFunction(Q)
du = TrialFunction(Q)
form = phi*du*dP
q = assemble(form)
I would like this example to produce the identity matrix. I understand that if dP does imply an integral for this to be mathematically correct, the basis functions would have to be delta functions, but is there a way to just consider this as a sum rather than an integral over a point? It seems that being able to use dP in this way has been discussed in a few blueprints, and it also seems to behave correctly when used in the context of the PointIntegral multistage solver.
Hand coding a point source by manipulating vectors is fairly straightforward. My problem has a time derivative on the point terms, so I need to be able to use dP in conjunction with test and trial functions. Is there a way to do this?
