dolfin team mailing list archive

[Anders Logg <logg@xxxxxxxxx>]

On Wed, Sep 17, 2008 at 01:23:17PM +0200, Evan Lezar wrote:
> 
> 
> On Tue, Sep 16, 2008 at 2:18 PM, Anders Logg <logg@xxxxxxxxx> 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()

-- 
Anders

Attachment: signature.asc
Description: Digital signature