← Back to team overview

dolfin team mailing list archive

Re: Function and DofMap

 

On Sun, Sep 07, 2008 at 09:15:40AM +0200, Andy Ray Terrel wrote:
> On Sun, Sep 7, 2008 at 12:41 AM, Martin Sandve Alnæs <martinal@xxxxxxxxx> wrote:
> > 2008/9/6 Anders Logg <logg@xxxxxxxxx>:
> >> On Sat, Sep 06, 2008 at 04:22:09PM +0200, Martin Sandve Alnæs wrote:
> >>> 2008/9/6 Garth N. Wells <gnw20@xxxxxxxxx>:
> >>> >
> >>> >
> >>> > Dag Lindbo wrote:
> >>> >> Anders Logg wrote:
> >>> >>> There seems to be a problem (among many) with the current design of
> >>> >>> the Function classes (see thread "evaluating higher order mesh function").
> >>> >>>
> >>> >>> In particular, the finite element is missing in DiscreteFunction. My
> >>> >>> suggestion would be to just add it and let a DiscreteFunction consist
> >>> >>> of the following four items which are always available:
> >>> >>>
> >>> >>>   mesh, x, dof_map, finite_element
> >>> >>>
> >>> >>> Is this enough, and what other issues to we need to fix?
> >>> >>>
> >>> >>
> >>> >> One major issue which I just want to reiterate is ownership of data. As
> >>> >> it stands, the DiscreteFunction may or may not be responsible for e.g.
> >>> >> the dof vector x, depending on whether local_vector is a NULL pointer or
> >>> >> not. Take a look at the thread "Ownership" from Garth on 06/26/2008.
> >>> >>
> >>> >
> >>> > Yes, this is a big problem and has caused me a few headaches with bugs.
> >>> > For example, passing a user-defined Function to a function to convert it
> >>> > to a DiscreteFunction via a projection onto a finite element basis
> >>> > causes a problem because the FiniteElement which the projected Function
> >>> > points to goes out of scope once the function is exited.
> >>> >
> >>> >> A problem related to this is initialization of the DiscreteFunction. We
> >>> >> had a bug previously where the LinearPDE class maintained ownership of
> >>> >> the solution vector. The only way to prevent this was to break the
> >>> >> encapsulation of DiscreteFunction by making it a friend of LinearPDE (as
> >>> >> with XMLFile for the same reasons). Here is some of the code that
> >>> >> handles this initializaton today (L101 in LinearPDE.cpp):
> >>> >>
> >>> >>   u.init(mesh, *x, a, 1);
> >>> >>   DiscreteFunction& uu = dynamic_cast<DiscreteFunction&>(*u.f);
> >>> >>   uu.local_vector = x;
> >>> >>
> >>> >> This ain't poetry in my opinion :)
> >>> >>
> >>> >
> >>> > Indeed, this isn't nice, and there is something similar in XMLFile.cpp.
> >>> >
> >>> > Garth
> >>>
> >>> We should start to use std::tr1::shared_ptr. There is some support for it
> >>> with python in swig 1.3.35, which is part of the upcoming Ubuntu Intrepid
> >>
> >> The main issue is how we want to initialize Functions, and if one
> >> should allow to set members.
> >>
> >> For simplicity, say that a Function is defined only by a Vector.
> >> Then we have a few different situations to consider:
> >>
> >> 1. Function creates the Vector
> >>
> >>   Function u;
> >>   Vector& x = u.vector();
> >>
> >> 2. Function gets the Vector
> >>
> >>   Vector x;
> >>   Function u(x);
> >>
> >> 3. Function gets initialized with a Vector
> >>
> >>   Function u;
> >>   Vector x;
> >>   u.init(x);
> >>
> >> Do we want to support all of 1-3? Things become considerable easier if
> >> we can make some simplifying assumptions.
> >>
> >> How visible would a shared_ptr be in the interface?
> >
> > A shared_ptr must be visible to the user every single place
> > a pointer is passed around, otherwise the reference count
> > won't be correct and we'll just have more problems.
> >
> > Having different options for ownership of data isn't just
> > an interface issue, it's important for coupling dolfin code
> > with other code.
> >
> > It's also not only Function that's affected by this issue,
> > you might have a PETSc Vec or Epetra_FEVector and
> > wish to assemble into if, so you'll need to wrap it in a
> > dolfin::*Vector during assembly, without giving up
> > ownership of the data.
> >
> >
> 
> As far as an interface is concerned I like having a vector hidden from
> me so only having the 2nd option is fine, if all the other
> functionality is working as it should.  It might just be because of
> all the functional programming UoC has made me teach but I really
> don't like the init methods in the user interface.  My preference is
> to have the copy constructors and assignment operators work correctly
> so that to reinitialize one just needs to do something like:
> 
> Function u;
> Vector x;
> u = Function(x);

I think that looks good, but the last line in your example will first
call the constructor in Function and then everything will be copied
using the assignment operator. But I guess this is not a problem if
the members are shared pointers.

I also dislike the init() functions and I try to avoid them, but very
often we need to add them in to do things like reading from file:

  Function u:
  File file("solution.xml");
  file >> u; // init() function called in XML reader

> BTW, it seems this discussion on function and dof maps has been going
> on for a while and spread out over many email threads.  Would this be
> a good issue to keep track of on bugzilla?
> 
> -- Andy

Would that work? I feel it's easier to have a discussion on the
mailing list than on the bug tracker. The problem is that all the
times we've been discussing this, we've lost momentum after some time
and then we forget about it. Let's see if we can finish up this time.

-- 
Anders

Attachment: signature.asc
Description: Digital signature


Follow ups

References