Martin Sandve Alnæs wrote:
> 2008/9/8 Garth N. Wells <[EMAIL PROTECTED]>:
>>
>> Anders Logg wrote:
>>> On Sun, Sep 07, 2008 at 03:38:44PM +0100, Garth N. Wells wrote:
>>>> Anders Logg wrote:
>>>>> 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 <[EMAIL PROTECTED]>:
>>>>>>>> On Sun, Sep 07, 2008 at 08:27:41AM +0100, Garth N. Wells wrote:
>>>>>>>>> Martin Sandve Alnæs wrote:
>>>>>>>>>> 2008/9/6 Anders Logg <[EMAIL PROTECTED]>:
>>>>>>>>>>> On Sat, Sep 06, 2008 at 04:22:09PM +0200, Martin Sandve Alnæs wrote:
>>>>>>>>>>>> 2008/9/6 Garth N. Wells <[EMAIL PROTECTED]>:
>>>>>>>>>>>>> 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.
>>>>>
>>>> I think that having two options would lead to confusion.
>>>>
>>>> On a another note, I don't really like the above code anyway. We could
>>>> avoid this in many cases by adding the constructor
>>>>
>>>> Function::Function(Mesh& mesh, Form& form, uint i = 1);
>>>>
>>>> and letting Function create the Vector. Function already has a member
>>>> function for accessing the underlying Vector.
>>>>
>>>> Garth
>>> That's what I suggested earlier. It would simplify many things if we
>>> always know who owns what. We could assume that the Function never
>>> owns the Mesh but always owns the Vector (for discrete functions).
>>>
>>> If someone has a Vector and wants to create a Function from it, then
>>> just do
>>>
>>> Vector x;
>>> Function u;
>>> u.vector() = x; // Copy data
>>>
>> OK, and if someone really wants a Function to point to x, they can do
>>
>> shared_ptr<Vector> x(new Vector);
>> Function u;
>> u.vector() = x; // Share data
>>
>> If Function uses shared_ptr internally for its Vector, it doesn't have
>> to worry about keeping track of ownership.
>>
>> Garth
>
> You can't overload a function on its return type.
>
> For this to work:
>>> u.vector() = x; // Copy data
> vector() must return a reference to the Vector.
>
> For this to work:
>> u.vector() = x; // Share data
> vector() must return a reference to the shared_ptr.
>
> What's wrong with just using a set-function?
>
> shared_ptr<Vector> v = u.vector();
>
> shared_ptr<Vector> shared_ptr_x;
> u.setVector(shared_ptr_x); // Share data
>
> Vector vector_x;
> u.setVector(vector_x); // Copy data
>
Looks ok. Only issue is aesthetic - we try to avoid 'set' and 'get'
functions.
Garth
>
> Also, if u.vector returns a shared_ptr, you can't assign it to a
> reference, but you should be able to dereference it:
>
> Vector & v = *u.vector();
>
> --
> Martin
_______________________________________________
DOLFIN-dev mailing list
[email protected]
http://www.fenics.org/mailman/listinfo/dolfin-dev
