dolfin team mailing list archive

[Anders Logg <logg@xxxxxxxxx>]

On Wed, Oct 22, 2008 at 07:45:03PM +0100, Garth N. Wells wrote:
> 
> 
> Anders Logg wrote:
> > On Wed, Oct 22, 2008 at 08:24:18PM +0200, Anders Logg wrote:
> >> On Wed, Oct 22, 2008 at 07:05:09PM +0100, Garth N. Wells wrote:
> >>>
> >>> Anders Logg wrote:
> >>>> On Wed, Oct 22, 2008 at 12:08:43PM +0100, Garth N. Wells wrote:
> >>>>> The thread on the new Function design has digressed from the immediate 
> >>>>> issue, so I'm restarting it.
> >>>>>
> >>>>> The issue is how to deal with user defined functions in C++. What if we 
> >>>>> have a design such that:
> >>>>>
> >>>>> - All Functions must have a FunctionSpace
> >>>> Yes.
> >>>>
> >>>>> - A FunctionSpace does not have to be complete (a complete FunctionSpace 
> >>>>> having a Mesh, a FiniteElement and a DofMap). As a minimum requirement, 
> >>>>> and FunctionSpace must have a Mesh.
> >>>> I think it should always need to be complete. We can avoid a lot of
> >>>> trouble if we make firm requirements: a Function is always associated
> >>>> with a FunctionSpace and the FunctionSpace is always completely
> >>>> defined.
> >>>>
> >>>>> - A FiniteElement and/or DofMap can be attached to a FunctionSpace after 
> >>>>> its creation
> >>>>>
> >>>>>
> >>>>> Related to the functions Function::interpolate(double* coefficients, ..) 
> >>>>> for interpolating Functions on cells
> >>>>>
> >>>>> - The new Function::interpolate functions do not take a FiniteElement as 
> >>>>> an argument, so it is not possible to interpolate a function in a 
> >>>>> different space. Is it desirable to allow Functions from one space to be
> >>>>> interpolated in another?
> >>>> Yes, it's desirable. This can be done globally when needed:
> >>>>
> >>>>   Function v(V);
> >>>>   Vector coefficients;
> >>>>   v.interpolate(coefficients, W);
> >>>>
> >>>> This will compute the global coefficient vector for v on W.
> >>>>
> >>>> This is currently implemented with the assumption that the meshes for
> >>>> V and W are the same but could quite easily be extended to
> >>>> non-matching meshes now that functions can be evaluated at arbitrary
> >>>> points (using GTS).
> >>>>
> >>>>> If we do this, would:
> >>>>>
> >>>>> - The above allow FiniteElement types to be checked at runtime for 
> >>>>> consistency (the FiniteElement passed to Function::interpolate should be 
> >>>>> the same as the Functions own FiniteElement for discrete Functions. This 
> >>>>> is what we did with the old design.)>
> >>>>>
> >>>>> - The above deal with the issue of user-defined functions which have a a 
> >>>>> FunctionSpace but no FiniteElement?
> >>>>>
> >>>>> - The above deal with special functions, like the mesh size h?
> >>>> Here's my suggestion for how to handle initialization of Functions in C++.
> >>>> There is no need for a circular dependency. First some simple facts:
> >>>>
> >>>>   1. The constructor of a Form may require one or more Functions
> >>>>   2. The constructor of a Function requires a FunctionSpace
> >>>>
> >>>> From this it follows that we must do something like
> >>>>
> >>>>   Function f(V);
> >>>>   Poisson::BilinearForm a;
> >>>>   Poisson::LinearForm L(f);
> >>>>
> >>>> The question is now how to initialize the FunctionSpace V. My
> >>>> suggestion would be to extend the code generation to generate code for
> >>>> creating the FunctionSpace(s), just like we do already when we
> >>>> generate UFC code + some extra code for defining the Form classes
> >>>> (when using -l dolfin in FFC). This does not make FFC DOLFIN-specific
> >>>> or DOLFIN FFC-specific. One can still use other form compilers, but the
> >>>> interface will not be as nice.
> >>>>
> >>>> So, one would do something like this:
> >>>>
> >>>>   #include "Poisson.h"
> >>>>
> >>>>   int main()
> >>>>   {
> >>>>     Mesh mesh("mesh.xml");
> >>>>     Poisson::TestSpace V(mesh);
> >>> Why 'TestSpace'?
> >>>
> >>> Also, when many functions are present, how would we 
> >>> identify each one? Could FFC create a class Poisson::FunctionSpace::foo 
> >>> when 'foo' is the name of the function in the FFC input?
> >> We could let FFC create a number of function space classes:
> >>
> >>   Poisson::TestSpace
> >>   Poisson::TrialSpace
> >>   Poisson::FunctionSpace_0
> >>   Poisson::FunctionSpace_1
> >>   ...
> >>
> >> or even
> >>
> >>   Poisson::v::FunctionSpace
> >>   Poisson::u::FunctionSpace
> >>   Poisson::f::FunctionSpace
> >>   Poisson::g::FunctionSpace
> >>
> >> (but I'm not sure it looks very nice).
> >>
> >> A complication I didn't think of before is how to handle the case when
> >> all spaces are the same (to get reuse of dofmaps). One option would be
> >> to either generate just one class if all spaces are the same, and
> >> otherwise generate separate clases like above. So either one does
> >>
> >>   Poisson::FunctionSpace V(mesh);
> >>   Function f(V);
> >>   Function g(V);
> >>   Poisson::LinearForm(f, g);
> >>
> >> or
> >>
> >>   Poisson::FunctionSpace_0 V(mesh);
> >>   Poisson::FunctionSpace_1 W(mesh);
> >>   Function f(V);
> >>   Function g(W);
> >>   Poisson::LinearForm(f, g);
> >>
> >>>>     Function f(V);
> >>>>     Poisson::BilinearForm a;
> >>>>     Poisson::LinearForm L(f);
> >>>>     ...
> >>>>   }
> >>>>
> >>>> First V, then f, then L.
> >>>>
> >>> Would this work for all the Functions in SpecialFunction.h?
> >> I haven't thought about it, but yes I guess it would.
> >>
> >> We can provide a set of predefined FunctionSpaces (like DG on
> >> triangles and tetrahedra).
> > 
> > I'm moving the discussion of the removal of DofMapSet over here since
> > it's related.
> > 
> > If a Form should own a set of FunctionSpaces (replacing the current
> > DofMapSet), that would prevent reuse of FunctionSpaces (and DofMaps)
> > across multiple Forms.
> > 
> > Another option would be to always require that a Form is initialized
> > with one or more FunctionSpaces. This makes sense and is in agreement
> > with the requirement of a FunctionSpace when initializing a Function.
> > 
> > It would then be
> > 
> >   PoissonFunctionSpace V(mesh);
> >   Function f(V);
> >   Function g(V);
> >   PoissonBilinearForm a(V, V);
> >   PoissonLinearForm L(V, f, g);
> > 
> > Then it's possible to reuse V across multiple forms (even forms not
> > included in Poisson.h).
> > 
> 
> Could work. To keep things clear cut, I would use
> 
>     PoissonFunctionSpace V(mesh);
>     Array<FunctionSpace*> VxV(V, V);
>     PoissonBilinearForm a(VxV);  // Might also nee a shared_ptr version
>     PoissonLinearForm L(V, f, g);
> 
> The first argument to a Form defines the test and/or trial spaces, the 
> other arguments are the Functions in the Form.
> 
> Garth

That looks like a complication. The constructor of a Form takes a
variable number of arguments anyway, so it's no big deal to require
one FunctionSpace for a linear form and two for a bilinear form.
It also looks nicer (avoiding Array and pointers).

-- 
Anders

Attachment: signature.asc
Description: Digital signature