dolfin team mailing list archive

[Anders Logg <logg@xxxxxxxxx>]

On Mon, Jul 07, 2008 at 10:38:52AM +0200, Martin Sandve Alnæs wrote:
> 2008/7/7 Anders Logg <logg@xxxxxxxxx>:
> > On Fri, Jul 04, 2008 at 11:25:52AM +0200, Martin Sandve Alnæs wrote:
> >> 2008/7/3 Martin Sandve Alnæs <martinal@xxxxxxxxx>:
> >> > 2008/7/3 Martin Sandve Alnæs <martinal@xxxxxxxxx>:
> >> >> I've imple
> >> >>
> >> >> 2008/7/3 Dag Lindbo <dag@xxxxxxxxxx>:
> >> >>> Martin Sandve Alnæs wrote:
> >> >>>> 2008/7/2 Dag Lindbo <dag@xxxxxxxxxx>:
> >> >>>>> Martin Sandve Alnæs wrote:
> >> >>>>>> I just added checks for valid Functions to the assembler, but this
> >> >>>>>> turned out to break some demos so I changed it to warnings.
> >> >>>>> This seems like helpful checks to make!
> >> >>>>>
> >> >>>>>> The reason is that Function(mesh, 0.0), which creates a ConstantFunction,
> >> >>>>>> is used for zero vector functions as well, which leaves the value shape
> >> >>>>>> of such a function undefined. This makes error checking impossible,
> >> >>>>>> which is not very nice.
> >> >>>>>>
> >> >>>>>> In my opinion, the correct next step is to keep these warnings,
> >> >>>>>> and turn them into errors after "a while". We should either add
> >> >>>>>> support for ConstantFunctions with tensor shape (would be nice),
> >> >>>>>> or make the demos and apps use user-define functions.
> >> >>>>> Adding tensor shape to the existing concept of ConstantFunction seems
> >> >>>>> like the way to go. IMHO, the convenience if this over user defined
> >> >>>>> functions is too great to discard.
> >> >>>>>
> >> >>>>> Would it not just amount to having ConstantFunction::rank() return the
> >> >>>>> proper value as determined by the mutable member size, instead of 0? (on
> >> >>>>> line 28 of ConstantFunction.cpp)
> >> >>>>
> >> >>>> No, because this size isn't known at construction time, only when
> >> >>>> interpolate is called with a finite_element object. A suitable (set
> >> >>>> of) constructor(s) needs to be defined, mirrored in Function, and
> >> >>>> handled in the swig interface for python.
> >> >>>
> >> >>> Yes, I see that this is needed if the rank and shape of the function is
> >> >>> to be set constructor. Clearly, this makes for a nice user interface! My
> >> >>> suggestion was more along the lines of a quick fix to get rid of the
> >> >>> warnings: rank() returns zero if size == 1, returns 1 if size == 2 or 3 etc.
> >> >>>
> >> >>> /Dag
> >> >>>
> >> >>>>
> >> >>>> // Scalar constructor
> >> >>>> Function f(mesh, 0.0);
> >> >>>>
> >> >>>> // General tensor constructor:
> >> >>>> Function f(Mesh& mesh, uint rank, uint* shape, uint* values);
> >> >>>>
> >> >>>> // vector
> >> >>>> uint rank = 1;
> >> >>>> uint shape[1] = { 2 };
> >> >>>> real values[2] = { 1., 2. };
> >> >>>> Function f(mesh, rank, shape, values);
> >> >>>> // or
> >> >>>> uint size = 2;
> >> >>>> real values[2] = { 1., 2. };
> >> >>>> Function f(mesh, 1, &size, values);
> >> >>>>
> >> >>>> // 3x3 tensor:
> >> >>>> uint rank = 2;
> >> >>>> uint shape[2] = { 3, 3 };
> >> >>>> real values[9] = {
> >> >>>>    1., 2., 3.,
> >> >>>>    4., 5., 6.,
> >> >>>>    7., 8., 9. };
> >> >>>> Function f(mesh, rank, shape, values);
> >> >>>>
> >> >>>
> >> >>>
> >> >>
> >> >> I just implemented this, but can't reproduce the old behaviour at the
> >> >> same time, so the demos must be updated and I guess many user
> >> >> applications as well. I'd like to have consensus to do this change
> >> >> before I continue with the demos and check in.
> >> >>
> >> >> The typical example is a zero function intended to be a vector:
> >> >>  Function f(mesh, 0.0);
> >> >>
> >> >> Which now would become:
> >> >>  uint fsize = 3;
> >> >>  real fvalues[3] = { 0.0, 0.0, 0.0 };
> >> >>  Function f(mesh, 1, &fsize, fvalues);
> >> >>
> >> >> I could evt. add a vector constructor:
> >> >>  real fvalues[3] = { 0.0, 0.0, 0.0 };
> >> >>  Function f(mesh, 3, fvalues);
> >> >
> >> > And:
> >> >
> >> >  Function f(mesh, 3, 0.0);
> >> >
> >> > which does the same as the above.
> >> >
> >>
> >> I'd rather not wait a month before I merge, so I'll push these changes.
> >> Somebody who knows the logic should look at dolfin_function_pre.i
> >
> > I think Ola knows these typemaps best.
> >
> >> and fix the typemaps there, until then the python syntax for constant
> >> vector functions will be:
> >>
> >> f = Function(mesh, size, array)
> >>
> >> E.g.
> >> f = Function(mesh, 3, numpy.array((1.,2.,3.)))
> >>
> >> Or simply:
> >> f = Function(mesh, 3, 0.0)
> >
> > I think that it would be nicer to create a vector-valued function in
> > C++ by
> >
> >  Array values(0.0, 0.0, 0.0);
> >  Function f(mesh, values);
> 
> Ok! I wasn't sure which array class we're supposed to use in C++ dolfin...

The reason we have our own Array class (which is just std::vector) is
to avoid confusion with two classes vector and Vector.

This is mostly important in interfaces. We use std::vector and other
STL classes freely internally.

-- 
Anders

Attachment: signature.asc
Description: Digital signature