dolfin team mailing list archive

["Martin Sandve Alnæs" <martinal@xxxxxxxxx>]

2008/11/3 Anders Logg <logg@xxxxxxxxx>:
> On Mon, Nov 03, 2008 at 10:36:52AM +0100, Martin Sandve Alnæs wrote:
>> 2008/11/3 Johan Hake <hake@xxxxxxxxx>:
>> > On Monday 03 November 2008 09:00:58 Martin Sandve Alnæs wrote:
>> >> 2008/11/2 Johan Hake <hake@xxxxxxxxx>:
>> >> > On Sunday 02 November 2008 17:44:39 Anders Logg wrote:
>> >> >> On Sun, Nov 02, 2008 at 04:28:45PM +0100, Johan Hake wrote:
>> >> >> > Hello!
>> >> >> >
>> >> >> > I have fixed the swig interface so it compiles. I had to add
>> >> >> > SubFunction in dolfin_function.h. If we do not want this we have to
>> >> >> > add
>> >> >> >
>> >> >> >   #include "SubFunction.h"
>> >> >> >
>> >> >> > to Function.h, as it is in 0.8.1.
>> >> >> >
>> >> >> > As function inherits ufc::function I also %imported ufc.h and added
>> >> >> > ufc-1 as swig dependencies in scons.cfg.
>> >> >> >
>> >> >> > I %renamed function.in to function._in and %extended the interface to
>> >> >> > SubFunction so one can write
>> >> >> >
>> >> >> >   u in V
>> >> >> >
>> >> >> > as discussed previously.
>> >> >>
>> >> >> Very nice!
>> >> >>
>> >> >> > I am curious of how we will fix the whole PyDOLFIN interface with the
>> >> >> > precompiled function spaces aso. These can all be hidden in the
>> >> >> > assemble function and/or in the LinearPDE class. But what if we want
>> >> >> > to expose a FunctionSpace from a compiled form? Do we want this?
>> >> >> >
>> >> >> > Or what if we want to create a FunctionSpace from an element and a
>> >> >> > mesh. This can be done more or less in the same way as we have solved
>> >> >> > Function in 0.8.1, but then we miss the whole point of reusing
>> >> >> > FunctionSpaces.
>> >> >> >
>> >> >> > Any good suggestions?
>> >> >> >
>> >> >> > Johan
>> >> >>
>> >> >> Here's a good suggestion... :-)
>> >> >>
>> >> >> We create a new Python class (in site-packages/dolfin/functionspace.py)
>> >> >> which either inherits from FormCompiler.FiniteElement and
>> >> >> dolfin::FunctionSpace or owns such objects as member data.
>> >> >
>> >> > I suppose it have to inherit from FormCompiler.FiniteElement, but it only
>> >> > have to owe a cpp:FunctionSpace.
>> >>
>> >> To clarify: The (only?) reason for owning and not inheriting is that
>> >> we want to reuse function spaces.
>> >>
>> >> > So the __init__ function in FunctionSpace would be something like:
>> >> >
>> >> >  def __init__(self, mesh, family, order):
>> >> >      FormCompiler.FiniteElement.__init__\
>> >> >         (self,family,dim2shape(mesh.geometry().dim(),order))
>> >> >      if self.value_dimension(0) > 1:
>> >> >          form = ffc.TestFunction(self)[0]*ffc.dx
>> >> >      else:
>> >> >          form = ffc.TestFunction(self)*ffc.dx
>> >> >      (compiled_form, module, form_data) = jit(form,mesh)
>> >> >      self._V = compiled_form.function_space(0)
>> >> >
>> >> > and then we can add something like
>> >> >
>> >> >  def cpp_function_space(self):
>> >> >      return self._V
>> >> >
>> >> > to the interface, where the latter function is used in the __init__
>> >> > function of the Function class to initialize the cpp_Function?
>> >> >
>> >> > Maybee we for compatability should let FunctionSpace inherit
>> >> > cpp_FunctionSpace. If we do that the __init__ function from above would
>> >> > look a bit differently.
>> >>
>> >> To still allow reusing FunctionSpaces, it would have to be some kind
>> >> of envelope-letter design then.
>> >
>> > Can't it just pass it self to where ever it is needed, instead of passing the
>> > function space that it owe?
>>
>> That's the point of inheritance, yes. And that's why "itself" must act
>> as an "envelope" around the "letter" that it owns by forwarding calls,
>> giving an envelope-letter design.
>>
>> >> >> In PyDOLFIN, FunctionSpace will then be something that replaces both
>> >> >> FormCompiler.FiniteElement and dolfin::FunctionSpace:
>> >> >>
>> >> >>   V = FunctionSpace(mesh, "Lagrange", 1)
>> >> >>   v = TestFunction(V)
>> >> >>   u = TrialFunction(V)
>> >> >>   f = Function(V)
>> >> >>   a = dot(grad(v), grad(u))*dx
>> >> >>   L = v*f*dx
>> >> >>   A = assemble(a)
>> >> >>   b = assemble(L)
>> >> >>
>> >> >> Note that the "triangle" or "tetrahedron" argument is not needed since
>> >> >> this can be deduced from the mesh.
>> >> >
>> >> > Looks nice. Would it then be possible to reuse the cpp_FunctionSpace that
>> >> > is initialized for V when the form is compiled using jit?
>> >>
>> >> The reuse must happen when V is constructed, so this has nothing to do
>> >> with compiling forms with jit.
>> >
>> > If V is constructed as suggested above, see __init__, jit needs to instantiate
>> > the functionspace from the compiled form module, and then use this to
>> > instantiate the form. (In this case where we only is interested in the
>> > actuall function space we do not need jit to instantiate the form though)
>>
>> That jit is used to compile a finite_element and dof_map in FunctionSpace
>> is a separate use of jit from compiling the forms a and L. The
>> ufl.Form instances
>> a and L won't see the dolfin.FunctionSpace, only the ufl.FiniteElement.
>> So the form is never instantiated with any function space.
>>
>> If we want to construct dolfin.Form objects in PyDOLFIN, we need to do
>> something more.
>>
>> > But when jit is called in the assemble function, it needs to instantiate the
>> > form using the functions spaces that are provided through the uncompiled form
>> > instead of instantiating them from the compiled form module.
>>
>> Yes. I think that can work out ok. When fetching ufl.FiniteElement objects
>> from a ufl.Form, those objects will actually be the dolfin.FunctionSpace objects
>> because of the inheritance above.
>>
>> formdata = ufl.FormData(a)
>> # formdata.elements == [V, ...]
>>
>> Thus the mesh will actually be implicitly defined through the function
>> spaces, and can be optional in assemble!
>>
>> So this should be possible:
>>
>> mesh = ...
>> V = FunctionSpace(mesh, "CG", 1)
>> f = Function(V)
>> a = f**2*dx
>> b = assemble(a)
>>
>> or this way for separately defined functions:
>>
>> mesh = ...
>> V = FunctionSpace(mesh, "CG", 1)
>> f = Function(V)
>> a = f**2*dx
>> fc = compiled_function...
>> b = assemble(a, {f: fc})
>>
>> However, I'd still like to be able to define forms separately in .ufl files
>> instead of always having to merge them with the PyDOLFIN application.
>> In that case we need to instantiate dolfin.FunctionSpace from the elements
>> in a given form. And we need to reuse those across forms, e.g. a and L.
>> (But at least then we can get the function names from the .ufl file.)
>>
>> forms = ufl.load_forms(ufl_filename) # can add returning functions later
>>
>> a, L = forms
>> formdata_a = ufl.FormData(a)
>> formdata_L = ufl.FormData(L)
>>
>> a_spaces = [dolfin.FunctionSpace(mesh, element) for element in
>> formdata_a.elements]
>> L_spaces = [dolfin.FunctionSpace(mesh, element) for element in
>> formdata_a.elements]
>>
>> a_form = dolfin.Form(a, a_spaces)
>>
>> Eh... This is getting complicated.
>> Maybe we can generate PyDOLFIN code as well ;)
>
> Why do you need to put the forms in separate .ufl files? It would be
> simpler to put them in separate .py files if you want them separate
> and use FunctionSpace (from PyDOLFIN).

For nontrivial equations it would be nice to be able to reuse the form
in both C++ and Python. It's a bit strange if we call UFL a language
and can't actually use it in Python. If the form is in a .ufl file it can
be analysed and debugged using commandline UFL utilities, and
once UFL stabilizes so will the implementation of the PDE.

Anyway, putting it in a .py file is close enough for the pragmatic
version of me, since it will only differ by some initial element definitions.


For pure UFC assembly this should work out fine:

ufc_a = ufc::form from somewhere

V = FunctionSpace(mesh, ...) # need to manually match this with ufc_a
f = Function(V, ...)

# either:
a = Form(ufc_a, V, V)
A = assemble(mesh, a, coefficients = [f,...])
# or:
A = assemble(mesh, ufc_a, spaces=[V, V], coefficients = [f,...])

If we add coefficient names to ufc::form we can allow named functions
here as well.

--
Martin