dolfin team mailing list archive

[Anders Logg <logg@xxxxxxxxx>]

On Sun, Sep 07, 2008 at 03:29:14PM +0200, Anders Logg wrote:
> On Sun, Sep 07, 2008 at 03:11:51PM +0200, Martin Sandve Alnæs wrote:
> > 2008/9/7 Anders Logg <logg@xxxxxxxxx>:
> > > On Sun, Sep 07, 2008 at 08:27:41AM +0100, Garth N. Wells wrote:
> > >>
> > >>
> > >> Martin Sandve Alnæs 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.
> > >> >
> > >>
> > >> So, in pseudo code, would it look something link this?
> > >>
> > >>    class DiscreteFunction
> > >>    {
> > >>      private:
> > >>
> > >>        shared_ptr<GenericVector> x;
> > >>
> > >>      public:
> > >>
> > >>        DiscreteFunction() : x(new Vector) {}
> > >>
> > >>        DiscreteFunction(shared_ptr<GenericVector> x)
> > >>        { x(x); }
> > >>
> > >>        shared_ptr<GenericVector> vec()
> > >>        {return x;}
> > >>    }
> > >> ?
> > >>
> > >> Garth
> > >
> > > What would the user code look like if we use shared_ptr for examples
> > > 1-3 above?
> > >
> > 
> > 
> > > >> 1. Function creates the Vector
> > > >>
> > > >>   Function u;
> > > >>   Vector& x = u.vector();
> > 
> > Function u;
> > Vector& x = u.vector(); // Storing this Vector& for later access is unsafe.
> > or
> > Function u;
> > shared_ptr<Vector> x = u.vector(); // Allows keeping the Vector around
> > after u is destroyed.
> > 
> > 
> > > >> 2. Function gets the Vector
> > > >>
> > > >>   Vector x;
> > > >>   Function u(x);
> > 
> > Vector x;
> > Function u(x); // Copy vector.
> >
> > shared_ptr<Vector> x = new Vector();
> > Function u(x); // Copy vector pointer, x or u may be deleted without
> > the other getting in trouble.
> 
> I don't think the first option is what one might expect, and I don't
> think the second example looks very nice.
> 
> We initialize Functions with a Mesh all the time and it would then be
> either very expensive to copy the mesh every time we create a Function
> from it (and one usually creates many functions on the same mesh), or
> we would have to write "shared_ptr" and "new" every time we used a
> Mesh.
> 
> Isn't there another option? I don't like the all the flags we have now
> like is_view, local_vector, etc, but this looks worse.

Is there a way to increase the count for a shared_ptr?

If there is, say a member named increase_ref(), we could do

  class DiscreteFunction
  {
  public:

      DiscreteFunction() : x(new Vector) {}

      DiscreteFunction(GenericVector& x) : x(x);
      {
          x.increase_ref(1);
      }

      DiscreteFunction(shared_ptr<GenericVector> x) : x(x) {}
      { x(x); }

      shared_ptr<GenericVector> vec()
      {return x;}

  private:

      shared_ptr<GenericVector> x;

  };

Then it would be possible to do

  Vector x;
  Function u(x);

as we can now and the Function u would know that someone else is
responsible for deleting the data.

Then if one writes code where the Vector goes out of scope, one must
use a shared_ptr, but not otherwise. We would not force everyone to
use shared_ptr all the time.

-- 
Anders

Attachment: signature.asc
Description: Digital signature