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
>
>
> ------------------------------------------------------------------------
>
> _______________________________________________
> DOLFIN-dev mailing list
> [email protected]
> http://www.fenics.org/mailman/listinfo/dolfin-dev
_______________________________________________
DOLFIN-dev mailing list
[email protected]
http://www.fenics.org/mailman/listinfo/dolfin-dev
