dolfin team mailing list archive
-
dolfin team
-
Mailing list archive
-
Message #09819
Re: PyDOLFIN Function
On Tue, Sep 23, 2008 at 05:07:49PM +0200, Evan Lezar wrote:
>
>
> On Tue, Sep 23, 2008 at 5:03 PM, Shawn Walker <walker@xxxxxxxxxxxxxxx> wrote:
>
>
>
> 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 <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()
> >>>
> >>>
> >>
> >>
> >> 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.
>
>
> I think that we are using this function to do different things. I am using the
> function to visualise a field based on a pre-computed vector of coefficients,
> whereas Anders is storing the solution of a set of linear equations directly in
> the function vector.
I hope things may clear up in a while, but we have many things going
right now: redesign of the function classes, the parallel XML parser
(and then the rest of the parallelization), introduction of the new
FunctionSpace class, higher-order meshes, renaming of member functions
etc. So it may be a few weeks until it clears up...
--
Anders
Attachment:
signature.asc
Description: Digital signature
References
