Add support for Dirac measure (point-source)

Registered by Anders Logg

The PointSource class can be replaced with a Dirac measure in the form compiler chain:

We set (formally)

  dp(\bar{x}) = \delta_{\bar{x}} dx

and then the generated code would need to

1. Check whether point is in the current triangle (easily handled with
   an appropriate code-snippet)

2. Evaluate the integrand at the point in question and add to the
   local tensor

This would require adding an addition integral to UFC:

  dirac_integral / point_integral / ???

Blueprint information

Series goal:
Milestone target:
Started by
Marie Rognes

Related branches


MSA: point_integral and *dP in UFC and UFL are in place in my personal branches.

One of the key design points is how to specify the actual point.

Currently the UFL syntax is u*dP(2), where 2 is a point ID analoguous to subdomain numbers.

Then DOLFIN somehow needs to associate this ID with an actual coordinate from e.g. a point set.
The advantage of this is that generated code is point- (and therefore mesh-) independent.

Actually, I think we should interpret u*dP(2) as the sum of point integrals over a point cloud. That will be useful for a set of measurement points, we wouldn't want to generate one integral class for each data point.

DHAM: Does it make more sense for that to be a sum or should it be regarded as a collection of integrals over a point cloud which therefore returns a set of values corresponding to the integral at each point?

Currently the UFC signature takes a cell (which DOLFIN must find), and a point x (passed from DOLFIN).

What's missing is FFC support and DOLFIN support. Volunteers? :)

AL: Point cloud sounds like a good solution. It is very compatible with the way we handle cell and facet integrals (in that it should be supplied from the DOLFIN side).

AL: A natural extension of the interface subdomain.mark_cells(mesh), subdomain.mark_facets() would be point.mark(mesh). The list of marked points could be stored in MeshDomains as part of the Mesh. Then no extra data will need to be passed to the assembler.

MSA: It might be a good idea to cache which cells these points reside in in the mesh data as well? But that is an optimization.

JH: Sounds cool. Stupid questions what happens when the point incidence with the vertex? Then that integral should only contribute to the degrees of freedom that also incidence with that vertex.

MER: I would like to have dP(x) where x is a Point or tuple in (Python) DOLFIN

MSA: dP(x) can be transformed in PyDOLFIN into a point cloud with x and dP(1) for the form compiler. Even dP([x1, x2, x3, x4, x5]) could be possible to do.

DHAM: I would suggest using dP[x] to signify a point. This makes a clean distinction between the syntax for point cloud IDs and for points themselves. It's also consistent with using dt[t] and dt[t1:t2] for a point in time and an interval in time in a future time measure (and indeed this is the notation which is currently proposed in the dolfin-adjoint time forms).

MER/PEF: What should be the default "domain" for dP in the assembler? dx assembles over all cells, dS/ds over all facets etc. Some options for dP would be: all vertices or all dof coordinates (for those elements defined via such)

MSA: All dof coordinates of which element? There is no element associated with the measure, and test and trial spaces can be different.

MSA: We may want to generate two integral functions in the point integral, one for arbitrary points in a given cell, and one for vertices/dof coordinates as mentioned above.

FFC can now generate code for point measures. How this is to be integrated in the assemblers is still up for discussion.

If dP defaults to integration over all vertices, it could immediately be integrated in
the assembler(s) in the same way as for the other integrals: by interating over
all vertices and adding contributions. This would then enable point sources
defined at vertices.

If dP is to follow some other regime, for instance involving a subset of arbitrary points, FFC needs to support point_integral::tabulate_tensor for arbitrary points (support for which awaits uflacs integration.)


Work Items

This blueprint contains Public information 
Everyone can see this information.


No subscribers.