On Tue, 23 Sep 2008, Anders Logg wrote:

> On Tue, Sep 23, 2008 at 12:35:48PM +0200, Marie Rognes wrote:
>> Anders Logg wrote:
>>> On Wed, Sep 17, 2008 at 01:23:17PM +0200, Evan Lezar wrote:
>>>
>>>> On Tue, Sep 16, 2008 at 2:18 PM, Anders Logg <[EMAIL PROTECTED]> wrote:
>>>>
>>>>     On Tue, Sep 16, 2008 at 01:10:35PM +0100, Garth N. Wells wrote:
>>>>    >
>>>>    >
>>>>    > Anders Logg wrote:
>>>>    >> On Mon, Sep 15, 2008 at 08:09:43PM +0100, Garth N. Wells wrote:
>>>>    >>>
>>>>    >>> Anders Logg wrote:
>>>>    >>>> On Mon, Sep 15, 2008 at 06:31:56PM +0200, Anders Logg wrote:
>>>>    >>>>> On Mon, Sep 15, 2008 at 04:42:32PM +0100, Garth N. Wells wrote:
>>>>    >>>>>> Anders Logg wrote:
>>>>    >>>>>>> On Mon, Sep 15, 2008 at 04:32:38PM +0100, Garth N. Wells wrote:
>>>>    >>>>>>>> Anders Logg wrote:
>>>>    >>>>>>>>> On Mon, Sep 15, 2008 at 03:47:59PM +0100, Garth N. Wells 
>>>> wrote:
>>>>    >>>>>>>>>> Anders Logg wrote:
>>>>    >>>>>>>>>>> On Mon, Sep 15, 2008 at 03:12:58PM +0100, Garth N. Wells
>>>>     wrote:
>>>>    >>>>>>>>>>>> Anders Logg wrote:
>>>>    >>>>>>>>>>>>> On Mon, Sep 15, 2008 at 03:06:55PM +0200, Johan Hake 
>>>> wrote:
>>>>    >>>>>>>>>>>>>> On Monday 15 September 2008 14:29:08 Garth N. Wells 
>>>> wrote:
>>>>    >>>>>>>>>>>>>>> Could a Python expert take a look at 
>>>> site-packges/dolfin/
>>>>     function.py?
>>>>    >>>>>>>>>>>>>>> The code directly following the comment
>>>>    >>>>>>>>>>>>>>>
>>>>    >>>>>>>>>>>>>>>   # Special case, Function(element, mesh, x), need to
>>>>     create simple form
>>>>    >>>>>>>>>>>>>>> to get arguments
>>>>    >>>>>>>>>>>>>>>
>>>>    >>>>>>>>>>>>>>> need to be updated but I don't understand it well.
>>>>    >>>>>>>>>>>>>> The first special case is for initializing a Function 
>>>> with
>>>>     a given Vector, by
>>>>    >>>>>>>>>>>>>> constructing a dofmap from the handed element.
>>>>    >>>>>>>>>>>>>>
>>>>    >>>>>>>>>>>>>> As constructing a Function from a vector is removed from
>>>>     the cpp interface,
>>>>    >>>>>>>>>>>>>> and we have not, (or have we?) figured out how to wrap a
>>>>     shared_ptr in swig,
>>>>    >>>>>>>>>>>>>> we should probably just remove the first case for now.
>>>>    >>>>>>>>>>>>>>
>>>>    >>>>>>>>>>>>>> Johan
>>>>    >>>>>>>>>>>>> The question is how we want to create discrete Functions 
>>>> in
>>>>     Python.
>>>>    >>>>>>>>>>>>> Previously, this was done by
>>>>    >>>>>>>>>>>>>
>>>>    >>>>>>>>>>>>>   u = Function(element, mesh, Vector())
>>>>    >>>>>>>>>>>>>
>>>>    >>>>>>>>>>>>> but now the third argument is not needed anymore. If we
>>>>     remove it,
>>>>    >>>>>>>>>>>>> we get
>>>>    >>>>>>>>>>>>>
>>>>    >>>>>>>>>>>>>   u = Function(element, mesh)
>>>>    >>>>>>>>>>>>>
>>>>    >>>>>>>>>>>>> but that doesn't work since that is the way to initialize 
>>>> a
>>>>    >>>>>>>>>>>>> user-defined function (something overloading eval()).
>>>>    >>>>>>>>>>>>>
>>>>    >>>>>>>>>>>>> We could put in a flag and make "discrete" the default. 
>>>> Then
>>>>     all
>>>>    >>>>>>>>>>>>> user-defined functions need to set the flag to "user".
>>>>    >>>>>>>>>>>>>
>>>>    >>>>>>>>>>>>> Suggestions? This is a good time to worry about how we 
>>>> want
>>>>     to design
>>>>    >>>>>>>>>>>>> the Function interface.
>>>>    >>>>>>>>>>>>>
>>>>    >>>>>>>>>>>> Sounds ok to me. This is basically what Vector() was doing,
>>>>     and a flag
>>>>    >>>>>>>>>>>> would be more descriptive.
>>>>    >>>>>>>>>>>>
>>>>    >>>>>>>>>>>> Garth
>>>>    >>>>>>>>>>> Maybe we could first try to think seriously about reducing 
>>>> the
>>>>     number
>>>>    >>>>>>>>>>> of different constructors in Function. There are 14 now! See
>>>>     below.
>>>>    >>>>>>>>>>>
>>>>    >>>>>>>>>>> I guess we need the following two basic constructors (empty
>>>>     and copy):
>>>>    >>>>>>>>>>>
>>>>    >>>>>>>>>>>   /// Create empty function (read data from file)
>>>>    >>>>>>>>>>>   Function();
>>>>    >>>>>>>>>>>
>>>>    >>>>>>>>>>>   /// Copy constructor
>>>>    >>>>>>>>>>>   Function(const Function& f);
>>>>    >>>>>>>>>>>
>>>>    >>>>>>>>>>> Then we have one for reading from file, which seems ok:
>>>>    >>>>>>>>>>>
>>>>    >>>>>>>>>>>   /// Create function from data file
>>>>    >>>>>>>>>>>   explicit Function(const std::string filename);
>>>>    >>>>>>>>>>>
>>>>    >>>>>>>>>>> And then the following set of constructors for constants:
>>>>    >>>>>>>>>>>
>>>>    >>>>>>>>>>>   /// Create constant scalar function from given value
>>>>    >>>>>>>>>>>   Function(Mesh& mesh, real value);
>>>>    >>>>>>>>>>>
>>>>    >>>>>>>>>> This one is useful.
>>>>    >>>>>>>>>>
>>>>    >>>>>>>>>>>   /// Create constant vector function from given size and
>>>>     value
>>>>    >>>>>>>>>>>   Function(Mesh& mesh, uint size, real value);
>>>>    >>>>>>>>>>>
>>>>    >>>>>>>>>> We could get rid of this one and use the below constructor.
>>>>    >>>>>>>>>>
>>>>    >>>>>>>>>>>   /// Create constant vector function from given size and
>>>>     values
>>>>    >>>>>>>>>>>   Function(Mesh& mesh, const Array<real>& values);
>>>>    >>>>>>>>>>>
>>>>    >>>>>>>>>> This one is useful.
>>>>    >>>>>>>>>>
>>>>    >>>>>>>>>>>   /// Create constant tensor function from given shape and
>>>>     values
>>>>    >>>>>>>>>>>   Function(Mesh& mesh, const Array<uint>& shape, const Array
>>>>     <real>& values);
>>>>    >>>>>>>>>>>
>>>>    >>>>>>>>>> This is the most generic of the constant functions, so I 
>>>> guess
>>>>     we need it.
>>>>    >>>>>>>>>>
>>>>    >>>>>>>>>>> And then there's this constructor which is needed for 
>>>> w.split
>>>>     (u, p):
>>>>    >>>>>>>>>>>
>>>>    >>>>>>>>>>>   /// Create discrete function from sub function
>>>>    >>>>>>>>>>>   explicit Function(SubFunction sub_function);
>>>>    >>>>>>>>>>>
>>>>    >>>>>>>>>>> But then there's the following mess of constructors:
>>>>    >>>>>>>>>>>
>>>>    >>>>>>>>>> Some of these constructors are necessary to support the
>>>>     PyDOLFIN
>>>>    >>>>>>>>>> interface. Can we get around this somehow to avoid 
>>>> duplication?
>>>>    >>>>>>>>>>
>>>>    >>>>>>>>>>>   /// Create function from given ufc::function
>>>>    >>>>>>>>>>>   Function(Mesh& mesh, const ufc::function& function, uint
>>>>     size);
>>>>    >>>>>>>>>>>
>>>>    >>>>>>>>>>>   /// Create discrete function for argument function i of 
>>>> form
>>>>    >>>>>>>>>>>   Function(Mesh& mesh, Form& form, uint i = 1);
>>>>    >>>>>>>>>>>
>>>>    >>>>>>>>>>>   /// Create discrete function for argument function i of 
>>>> form
>>>>    >>>>>>>>>>>   Function(Mesh& mesh, DofMap& dof_map, const ufc::form& 
>>>> form,
>>>>     uint i = 1);
>>>>    >>>>>>>>>>>
>>>>    >>>>>>>>>>>   /// Create discrete function for argument function i of 
>>>> form
>>>>     (data may be shared)
>>>>    >>>>>>>>>>>   Function(std::tr1::shared_ptr<Mesh> mesh,
>>>>    >>>>>>>>>>>            std::tr1::shared_ptr<GenericVector> x,
>>>>    >>>>>>>>>>>            std::tr1::shared_ptr<DofMap> dof_map, const
>>>>     ufc::form& form, uint i = 1);
>>>>    >>>>>>>>>>>
>>>>    >>>>>>>>>>>   /// Create discrete function based on signatures
>>>>    >>>>>>>>>>>   Function(std::tr1::shared_ptr<Mesh> mesh,
>>>>    >>>>>>>>>>>            const std::string finite_element_signature,
>>>>    >>>>>>>>>>>            const std::string dof_map_signature);
>>>>    >>>>>>>>>>>
>>>>    >>>>>>>>>>>   /// Create user-defined function (evaluation operator must
>>>>     be overloaded)
>>>>    >>>>>>>>>>>   explicit Function(Mesh& mesh);
>>>>    >>>>>>>>>>>
>>>>    >>>>>>>>>> We need this one.
>>>>    >>>>>>>>>>
>>>>    >>>>>>>>>> Garth
>>>>    >>>>>>>>> If we just consider discrete functions for a while, the 
>>>> question
>>>>     is
>>>>    >>>>>>>>> how these may be most conveniently (and naturally) defined in
>>>>     C++ and
>>>>    >>>>>>>>> Python.
>>>>    >>>>>>>>>
>>>>    >>>>>>>>> In C++, one only has a dolfin::Form, for example
>>>>     PoissonBilinearForm,
>>>>    >>>>>>>>> and then it's simple to create a discrete Function by
>>>>    >>>>>>>>>
>>>>    >>>>>>>>>   Function u(mesh, form);
>>>>    >>>>>>>>>
>>>>    >>>>>>>>> This will extract the element and dof map for the second
>>>>     argument of
>>>>    >>>>>>>>> the form (the trial function) which is normally what is 
>>>> needed.
>>>>    >>>>>>>>>
>>>>    >>>>>>>>> In Python, one does not have a dolfin::Form, but instead one 
>>>> has
>>>>     a
>>>>    >>>>>>>>> FiniteElement, and then the simplest thing to do is
>>>>    >>>>>>>>>
>>>>    >>>>>>>>>   u = Function(element, mesh)
>>>>    >>>>>>>>>
>>>>    >>>>>>>>> The element is the first argument for practical reasons (see
>>>>    >>>>>>>>> function.py) but maybe it shouldn't. I'd like to change this 
>>>> so
>>>>     that
>>>>    >>>>>>>>> the mesh is always first. All Functions require a Mesh and 
>>>> then
>>>>     it's
>>>>    >>>>>>>>> natural to put this first.
>>>>    >>>>>>>>>
>>>>    >>>>>>>>> So then we would have
>>>>    >>>>>>>>>
>>>>    >>>>>>>>>   C++:    Function u(mesh, form);
>>>>    >>>>>>>>>   Python: u = Function(mesh, element)
>>>>    >>>>>>>>>
>>>>    >>>>>>>>> On the other hand, we've been discussing adding a 
>>>> FunctionSpace
>>>>     class,
>>>>    >>>>>>>>> and then it might be natural to just have
>>>>    >>>>>>>>>
>>>>    >>>>>>>>>   C++:    Function u(V);
>>>>    >>>>>>>>>   Python: u = Function(V)
>>>>    >>>>>>>>>
>>>>    >>>>>>>>> This would create a discrete Function. Constant Functions and
>>>>    >>>>>>>>> user-defined Functions may be created without reference to a
>>>>    >>>>>>>>> FunctionSpace. This would solve the problem of overloading
>>>>    >>>>>>>>> constructors. It would be very clear that whenever a
>>>>     FunctionSpace is
>>>>    >>>>>>>>> involved, it is a discrete Function.
>>>>    >>>>>>>>>
>>>>    >>>>>>>> Agree. It's not appropriate to initialise a discrete function
>>>>     with a
>>>>    >>>>>>>> form. It seems that using a FunctionSpace will simplify the
>>>>     interface
>>>>    >>>>>>>> and provide uniformity across the C++ and Python interfaces, so
>>>>     let's
>>>>    >>>>>>>> get FunctionSpace (or something similar with another name) in
>>>>     place and
>>>>    >>>>>>>>   then remove some the Function constructors.
>>>>    >>>>>>> I think FunctionSpace is a good name.
>>>>    >>>>>>>
>>>>    >>>>>>>> Function spaces are not only associated with DiscreteFunctions.
>>>>     We
>>>>    >>>>>>>> usually interpolate user defined function in the finite element
>>>>     space,
>>>>    >>>>>>>> so perhaps there is some scope to unify discrete and user 
>>>> defined
>>>>    >>>>>>>> functions?
>>>>    >>>>>>>>
>>>>    >>>>>>>> Garth
>>>>    >>>>>>> Perhaps, but it would require storing an extra vector of values
>>>>     for
>>>>    >>>>>>> user-defined functions unnecessarily.
>>>>    >>>>>>>
>>>>    >>>>>> I'm not suggesting that we store a vector of values - just 
>>>> pointing
>>>>     out
>>>>    >>>>>> a function space is also associated with user-defined functions.
>>>>    >>>>>>
>>>>    >>>>>> Garth
>>>>    >>>>> Yes, I forgot for a minute that user-defined functions also need 
>>>> to
>>>>     be
>>>>    >>>>> associated with a particular function space when used in a form.
>>>>    >>>>>
>>>>    >>>>> Then it seems we still have the problem of differentiating between
>>>>    >>>>> user-defined functions and discrete functions when overloading
>>>>    >>>>> constructors.
>>>>    >>>> I have a suggestion that would solve this as well as another 
>>>> problem
>>>>    >>>> we've discussed a few times before regarding Functions
>>>>    >>>> (time-stepping). It's rather drastic but it might work.
>>>>    >>>>
>>>>    >>>> * Remove GenericFunction, DiscreteFunction, UserFunction etc and 
>>>> just
>>>>    >>>>   have one single class named Function.
>>>>    >>>>
>>>>    >>> Sounds good to me. I often found the current design a bit 
>>>> complicated.
>>>>     A
>>>>    >>> simpler design will make parallisation easier.
>>>>    >>>
>>>>    >>>> * Always require a FunctionSpace in the inialization of Function:
>>>>    >>>>
>>>>    >>>>   u = Function(V)
>>>>    >>>>
>>>>    >>>> * Let a Function consist of two things, a function space V and
>>>>    >>>>   degrees of freedom U:
>>>>    >>>>
>>>>    >>>>   u_h = \sum_i U_i \phi_i
>>>>    >>>>
>>>>    >>>>   where {\phi_i} are the basis functions for V.
>>>>    >>>>
>>>>    >>> As long as it's general enough to permit 'quadrature functions'.
>>>>    >>>
>>>>    >>>> * Keep a shared_ptr to V, and a pointer to a Vector U.
>>>>    >>>>
>>>>    >>> It might be useful to make the vector a shared pointer (need to 
>>>> check
>>>>    >>> what a smart pointer returns if it doesn't point to anything) 
>>>> because
>>>>     it
>>>>    >>> would be useful with sub-functions.
>>>>    >>>
>>>>    >>>> * The Vector pointer is initialized to 0 by default and as long as 
>>>> it
>>>>    >>>>   remains 0, the Function behaves like a user-defined Function, 
>>>> that
>>>>    >>>>   is, the eval() function is used.
>>>>    >>>>
>>>>    >>>> * Whenever someone calls the vector() member function, the Vector U
>>>>     is
>>>>    >>>>   initialized to a Vector of correct size (determined by the DofMap
>>>>     in
>>>>    >>>>   the FunctionSpace V). From this point on, the Function behaves 
>>>> like
>>>>    >>>>   a discrete Function, that is, the values in U are used, not 
>>>> eval().
>>>>    >>>>
>>>>    >>> Sounds good.
>>>>    >>>
>>>>    >>>> * This would make the behavior dynamic. For example, one may set u0
>>>>     to
>>>>    >>>>   an initial value in time-stepping and then do u0 = u1 and u0 
>>>> would
>>>>    >>>>   change behavior from user-defined to discrete.
>>>>    >>>>
>>>>    >>>> * Constant functions are handled by a new class named simply
>>>>     Constant,
>>>>    >>>>   which is a subclass of Function that overloads eval() and returns
>>>>    >>>>   the constant value.
>>>>    >>>>
>>>>    >>> Perhaps ConstantFunction would be a better name?
>>>>    >>
>>>>    >> Maybe, but it would also be natural to couple this to Constant in UFL
>>>>    >> in the same way as we couple Function to Function.
>>>>    >>
>>>>    >> Another reason to use Constant is that ConstantFunction implies that
>>>>    >> there may be other subclasses of Function than Constant and it would
>>>>    >> be good to avoid building a hierarchy (again). The class Constant 
>>>> just
>>>>    >> happens to be implemented as a subclass of Function, but it should be
>>>>    >> thought of as something different.
>>>>    >>
>>>>    >>> No matter what happens, it seems that FunctionSpace is the next 
>>>> step.
>>>>    >>> Related to this, is it the intention that a FiniteElement owns a
>>>>    >>> ufc::finite_elememt and provides wrapper, like DofMap?
>>>>    >>
>>>>    >> Yes, that was my plan. It only needs to overload the subset of
>>>>    >> functionality of ufc::finite_element that we need now.
>>>>    >>
>>>>    >
>>>>    > OK. I can put FiniteElement together pretty quickly if you haven't
>>>>    > started already.
>>>>    >
>>>>    > Garth
>>>>
>>>>     I haven't (not more than the empty template that's currently there).
>>>>
>>>>
>>>>
>>>> Hi
>>>>
>>>> What is the correct way to initialize a  discrete function in python.  I 
>>>> see
>>>> that the way I was using in the electromagnetic demo no longer works.  f =
>>>> Function(element, mesh, vector)
>>>>
>>>> Evan
>>>>
>>>
>>> That is the correct way, but the vector will be ignored since the
>>> Function from now on will create and own its vector of degrees of
>>> freedom.
>>>
>>> So you can do something like
>>>
>>>   f = Function(element, mesh, Vector())
>>>   x = f.vector()
>>>
>>> The vector you send in will be ignored and then you need to pick it
>>> out afterwards.
>>>
>>> Once the new design is in place it will be
>>>
>>>   f = Function(V) # V is a FunctionSpace
>>>   x = f.vector()
>>>
>>>
>>
>>
>> Follow-up question: So, what is the right way to initialize a Function
>> with a vector of computed values?
>
> You can't. You have to create the function first, then set the values
> in the vector:
>
>  f = Function(element, mesh, Vector())
>  x = f.vector()
>
>  solve(A, x, b)

I'm sorry, should this be:

f.vector() = x?

Maybe I just misread the question.

- Shawn
_______________________________________________
DOLFIN-dev mailing list
[email protected]
http://www.fenics.org/mailman/listinfo/dolfin-dev

Reply via email to