ffc team mailing list archive

[Kristian Oelgaard <k.b.oelgaard@xxxxxxxxxx>]

Quoting Johan Hake <hake@xxxxxxxxx>:

> On Wednesday 19 August 2009 16:18:48 Kristian Oelgaard wrote:
> > Quoting "Garth N. Wells" <gnw20@xxxxxxxxx>:
> > > Johan Hake wrote:
> > > > On Wednesday 19 August 2009 15:12:21 Garth N. Wells wrote:
> > > >> Johan Hake wrote:
> > > >>> On Wednesday 19 August 2009 14:57:20 Garth N. Wells wrote:
> > > >>>> Johan Hake wrote:
> > > >>>>> On Wednesday 19 August 2009 14:21:12 Garth N. Wells wrote:
> > > >>>>>> Garth N. Wells wrote:
> > > >>>>>>> Kent Andre wrote:
> > > >>>>>>>>> Not sure what we have to do actually. I suppose I want get any
> > > >>>>>>>>> clue from the
> > > >>>>>>>>>   demo/function/restriction
> > > >>>>>>>>>
> > > >>>>>>>>> demo?
> > > >>>>>>>>> I see that Kent has modified the dolfin.FunctionSpace to
> > > >>>>>>>>> reflect the present (broken) implementation of restriction,
> > > >>>>>>>>> anything in that line? Should such an implementation has been
> > > >>>>>>>>> done in FunctionSpaceBase instead of FunctionSpace? Kent?
> > > >>>>>>>>>
> > > >>>>>>>>> Johan
> > > >>>>>>>>
> > > >>>>>>>> This restriction is different in the sense that the restriction
> > > >>>>>>>> is made by the standard element on a sub mesh.
> > > >>>>>>>> The restriction in focus now is on the complete mesh but for
> > > >>>>>>>> only a part of the element.
> > > >>>>>>>> Combined, these two features will be very powerful.  I guess the
> > > >>>>>>>> sub-mesh restriction stuff is broken due to the parallel mesh.
> > > >>>>>>>
> > > >>>>>>> I've attached an example solver. To make it work for P^k with k >
> > > >>>>>>> 2 we need to restrict some functions to cell facets. Using FFC,
> > > >>>>>>> we do
> > > >>>>>>>
> > > >>>>>>>     scalar  = FiniteElement("Lagrange", "triangle", 3)
> > > >>>>>>>     dscalar = FiniteElement("Discontinuous Lagrange", "triangle",
> > > >>>>>>> 3)
> > > >>>>>>>
> > > >>>>>>>     # This syntax restricts the space 'scalar' to cell facets.
> > > >>>>>>>     mixed = scalar["facet"] + dscalar
> > > >>>>>>>
> > > >>>>>>> In PyDOLFIN we could do
> > > >>>>>>>
> > > >>>>>>>     V_dg = FunctionSpace(mesh, "DG", 3)
> > > >>>>>>>     V_cg = FunctionSpace(mesh, "CG", 3, "facet")
> > > >>>>>>>
> > > >>>>>>>     mixed = V_cg + V_dg
> > > >>>>>>
> > > >>>>>> I see now that we need to have the above option since the
> > > >>>>>> everything (e.g. the dof map) is created and compiled inside
> > > >>>>>> functionspace.py.
> > > >>>>>
> > > >>>>> We could implement it in FunctionSpaceBase, as an operator,
> > > >>>>> returning a new FunctionSpaceBase, with the restricted element as
> > > >>>>> arguemnt.
> > > >>>>>
> > > >>>>> However my concerns go to a sub space of a restricted mixed space
> > > >>>>> (if such exists) see the other email.
> > > >>>>>
> > > >>>>> If we go for 1, (as you have already done ;) ), we cannot create a
> > > >>>>> restricted mixed space using the addition syntax.
> > > >>>>>
> > > >>>>>   mixed = V1 + V2
> > > >>>>
> > > >>>> Do you mean to then restrict mixed? I'm already doing
> > > >>>
> > > >>> Precisely.
> > > >>>
> > > >>> After a second thought I think that doing it in the constructor is
> > > >>> more clean. Makes life simpler. Is this also possible i UFL, I mean:
> > > >>>
> > > >>>   e = FiniteElement("CG",cell,0,"facet")
> >
> > No, and the above element is not even legal, for CG the degree must be >=1
> > :)
> >
> > You have to use the __getitem__ operator of FiniteElementBase
> >
> > a = element[cell]
> >
> > or do
> >
> > e = ElementRestriction(element, domain)
> >
> > where domain is either a Measure or Cell and only Cell is supported in FFC.
> >
> > In FFC I might consider adding a default 'domain=None' to the argument list
> > of a FiniteElement since it will make some things easier. Maybe we should
> > add this in UFL too?
>
> More, harmony between the languages is nice, therefore I asked. Then we
> should
> probably support both syntaxes in PyDOLFIN, but the __getitem__ syntax might
> be troublesome to implement.

Personally I don't mind adding the 'domain' argument to FiniteElementBase in
UFL, but maybe someone have objections in this regard? Then at least this
syntax will be similar to what you will/can do in PyDOLFIN. The only difference
is then not being able to use __getitem__ in PyDOLFIN.

Kristian

> Johan
>
> > Kristian
> >
> > > >>> Should it?
> > > >>
> > > >> Is that a PyDOLFIN finite element? If it is, then yes, that's how it
> > > >> should work.
> > > >
> > > > I meant UFL.FiniteElement.
> > >
> > > OK. I don't know if that is the syntax. Kristian will know.
> > >
> > > >> Somewhat related, it's a bit annoying in PyDOLFIN that the code is
> > > >> recompiled when the mesh is changed. Can we avoid this? Should one
> > > >> then create a finite element and a dof map with which to initialise
> > > >> the Pythin FunctionSpace?
> > > >
> > > > Do you mean when one has refined a mesh?
> > >
> > > Yes, or just changed from one simulation to the next.
> > >
> > > As far as I understand it is the
> > >
> > > > reason for this that we need to recompute the dofmap.
> > >
> > > In C++ we create a dolfin::FunctionSpace from a dolfin::FiniteElement, a
> > > dolfin::DofMap and a dolfin::Mesh. The FiniteElement and DofMap are
> > > created with a ufc::finite_element and a ufc::dof_map, respectively.
> > >
> > > In the Python interface, in addition to what we have now, we would need
> > > FiniteElement and DofMap objects, with which we could create a
> > > FunctionSpace (together with a Mesh). Then, hopefully we could avoid
> > > recompilation when just the mesh is changed.
> > >
> > > Garth
> > >
> > > > How is this done in the C++ interface.
> > > >
> > > > Johan
> > > >
> > > >> Garth
> > > >>
> > > >>> Johan
> > > >>>
> > > >>>>      V0 = FunctionSpace(mesh, "CG", 3, "facet")
> > > >>>>      V1 = FunctionSpace(mesh, "DG", 3)
> > > >>>>
> > > >>>>      mixed = V0 + V1
> > > >>>>
> > > >>>> and it works. This type of restriction (we'll need to sort out our
> > > >>>> terminology) impacts the dofmap, which I why I think it needs to be
> > > >>>> involved in the creation of the function space.
> > > >>>>
> > > >>>> Garth
> > > >>>>
> > > >>>>> If we want a restricted mixed space we need to use the constructor.
> > > >>>>>
> > > >>>>>   mixed = MixedFunctionSpace((V1,V2),restriction="facet")
> > > >>>>>
> > > >>>>> Johan
> > > >>>>>
> > > >>>>>> Garth
> > > >>>>>>
> > > >>>>>>> or
> > > >>>>>>>
> > > >>>>>>>     V_dg = FunctionSpace(mesh, "DG", 3)
> > > >>>>>>>     V_cg = FunctionSpace(mesh, "CG", 3)
> > > >>>>>>>
> > > >>>>>>>     mixed = V_cg["scalar"] + V_dg
> > > >>>>>>>
> > > >>>>>>> Garth
> > > >>>>>>>
> > > >>>>>>>> Kent
> > > >>>>>>>>
> > > >>>>>>>> _______________________________________________
> > > >>>>>>>> DOLFIN-dev mailing list
> > > >>>>>>>> DOLFIN-dev@xxxxxxxxxx
> > > >>>>>>>> http://www.fenics.org/mailman/listinfo/dolfin-dev
> > > >>>>>>>
> > > >>>>>>> -----------------------------------------------------------------
> > > >>>>>>>---- -- -
> > > >>>>>>>
> > > >>>>>>> _______________________________________________
> > > >>>>>>> DOLFIN-dev mailing list
> > > >>>>>>> DOLFIN-dev@xxxxxxxxxx
> > > >>>>>>> http://www.fenics.org/mailman/listinfo/dolfin-dev
> > > >>>>>>
> > > >>>>>> _______________________________________________
> > > >>>>>> DOLFIN-dev mailing list
> > > >>>>>> DOLFIN-dev@xxxxxxxxxx
> > > >>>>>> http://www.fenics.org/mailman/listinfo/dolfin-dev
> > >
> > > _______________________________________________
> > > FFC-dev mailing list
> > > FFC-dev@xxxxxxxxxx
> > > http://www.fenics.org/mailman/listinfo/ffc-dev
> >
> > _______________________________________________
> > DOLFIN-dev mailing list
> > DOLFIN-dev@xxxxxxxxxx
> > http://www.fenics.org/mailman/listinfo/dolfin-dev
>