ffc team mailing list archive

[Anders Logg <logg@xxxxxxxxx>]

On Mon, May 23, 2005 at 01:48:55PM -1100, kirby@xxxxxxxxxxxx wrote:
> >I updated the TODO list with two suggested changes for FFC.
> >Tell me what you think.
> >
> >1. Projections
> >--------------
> >
> >I'm thinking of adding projections to the language. A 
> simple example:
> >
> >    P0 = FiniteElement("Discontinuous vector 
> Lagrange", "tetrahedron", 0)
> >    pi = Projection(P0)
> >
> >pi would then be a projection onto piecewise constants.
> >
> >This could be used to reduce the complexity of evaluating 
> certain forms.
> >If w is a Function, one can replace w with its projection 
> in a form:
> >
> >    a = pi(w)*u.dx(i)*v.dx(i)*dx
> >
> This seems like a good thing to do, and makes the 
> variational crime (projection) an explicit rather than 
> assumed thing in the language.  It is good to not have that 
> kind of semantic nuance.  Of course, we need a module for 
> generating quadrature loop nests in place if we want to 
> solve without the projection.  It would also be interesting 
> to have nodal interpolants and L2 projections as separate 
> things.  I don't know of any systematic explorations of 
> which are better in practice, so it might be worth doing.

A Function in FFC is always a member of a finite element space, and
each Function declared in FFC is always connected to a finite element
Function in DOLFIN.

The FFC Function is the restriction of a DOLFIN Function to a local
element.

Of course, one could change the syntax from

  f = Function(element)
  L = v*f*dx

to

  f  = Function()
  pi = Projection(element)
  L  = v*pi(f)*dx

but that would only make things more complicated to write (but perhaps
it would be more clear what happens).

What I meant was that in some cases, we want to project a Function
locally from one finite element space to another. Take as an example
the stabilization term of Navier-Stokes with linear elements on tets:

  a = U[j]*u[i].dx(j)*U[k]*v[i].dx(k)*dx

There are about 5 million integrals to be computed for the rank 8
reference tensor. The frozen velocity U is piecewise linear (obtained
from a previous solution). To speed up the compilation (and
evaluation) of the form, one can project U to piecewise constants,
which reduces the complexity by a factor 16. One can either do this
globally, by declaring a global piecewise constant function and doing
the projection in DOLFIN, or one can do it locally by telling FFC to
treat U as a piecewise constant function locally:

  pi = Projection(P0)
  a  = pi(U)[j]*u[i].dx(j)*pi(U)[k]*v[i].dx(k)*dx

Note that I have not implemented this yet.

> >2. Change .dx() --> .Dx()
> >-------------------------
> >
> >Ridg suggested changing the notation for derivatives 
> from .dx() to .Dx()
> >to avoid confusion with the integration operator *dx. What 
> do you think?
> >Poisson would then change from
> >
> >    a = u.dx(i)*v.dx(i)*dx
> >
> >to
> >
> >    a = u.Dx(i)*v.Dx(i)*dx
> >
> I agree that something like this is helpful as it avoids 
> overloading "dx".  But they differ only by case, so perhaps 
> we should consider:
> i.) u.deriv(i)
> or
> ii.) D(i)(u)

I think ".deriv()" is to long:

  u.deriv(i)*v.deriv(i)

And D(i) would look very much like a gradient:

  D(i)(u)*D(i)(v)

I need to think about this. At some point, we should add grad,
div, rot anyway.

> To be pedantic, in the second case, D is a higher order 
> function that takes the coordinate direction and returns the 
> function that differentiates a function in that direction.   

Yes.

/Anders