dolfin team mailing list archive

[Anders Logg <logg@xxxxxxxxx>]

On Sun, May 04, 2008 at 11:12:17PM +0200, Dag Lindbo wrote:
> 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

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