Subject: [FFC-dev] Assembly of singular source terms

From: "Dag Lindbo" <[email-masked]>
Date: Sun, 4 May 2008 23:12:17 +0200 (MEST)

Hello all!

I'm exploring how to use Dolfin and/or FFC to compute solutions to PDEs
involving singular source terms.

A typical model problem here is
-u_xx - u_yy = f*delta(gamma),
where f is some source density and delta is a Dirac distribution with
support in a set (curve) gamma in R^2. One has terms like this in
free-boundary flows etc. (The curve is unrelated to the eventual 2d mesh).

>From a mathematical point of view there is no problem integrating the
source term in every cell and assembling a contribution to the rhs in the
linear system. In fact, the element integral simplifies a line integral
that cuts through the cell (thus getting rid of the Dirac delta). How
could this be handled cleanly using FFC and Dolfin?

My (temporary?) solution: A secondary assembler in Dolfin. Code exists and
works (I'd be happy to share it on request).

The LineSourceAssembler class takes as input
  *) the curve as a set if ordered points (x_i,y_i)
  *) source density at each vertex, f_i
  *) the finite element from the linear form of the PDE
  *) mesh
and assembles the contribution into the rhs of the linear system. In a
nutshell, I do this by intersecting each line segment with the mesh (fast
using GTS), and manually computing the local element integrals in all
intersected cells. Some (small) amount of computation goes into computing
where and how the line segment cuts the cells. This only uses ffc (well
ufc I suppose) to evaluate basis functions. In my head I have a slightly
more elaborate design of a CustomAssembler abstract base class, where the
user provides a function that does the actual element integration.

However, I recognize that this is very un-Fenics. I would like to get some
input on how a better interface could be developed that ties it closer to
FFC and the core assembler.

Thanks!

/Dag

---

From: Anders Logg <[email-masked]>
Date: Sun, 4 May 2008 23:33:31 +0200

> FFC and the core assembler.

I don't have a clear view on how to design the interface to support
this type of integral, but one of my students (Kristoffer Selim) is
working on a similar problem where he needs to integrate terms (in his
case DG jump terms) over arbitrary curves cutting through the domain.
We have just started, but I hope in a few months time we'll know more.

The ideal thing would be to have something like

  void assemble(GenericTensor& A, Form& form, Mesh& curve);

-- 
Anders