dolfin team mailing list archive

[Anders Logg <logg@xxxxxxxxx>]

On Thu, Mar 12, 2009 at 05:41:42PM +0100, Martin Sandve Alnæs wrote:
> On Thu, Mar 12, 2009 at 5:26 PM, Anders Logg <logg@xxxxxxxxx> wrote:
> > On Thu, Mar 12, 2009 at 04:59:43PM +0100, Martin Sandve Alnæs wrote:
> >> On Thu, Mar 12, 2009 at 4:33 PM, Anders Logg <logg@xxxxxxxxx> wrote:
> >> > On Thu, Mar 12, 2009 at 01:42:55PM +0100, Martin Sandve Alnæs wrote:
> >> >> On Thu, Mar 12, 2009 at 1:16 PM, Johan Hake <hake@xxxxxxxxx> wrote:
> >> >> > On Thursday 12 March 2009 12:58:33 Garth N. Wells wrote:
> >> >> >> Anders Logg wrote:
> >> >> >> > On Thu, Mar 12, 2009 at 10:46:14AM +0100, Martin Sandve Alnæs wrote:
> >> >> >> >> On Thu, Mar 12, 2009 at 10:32 AM, Anders Logg <logg@xxxxxxxxx> wrote:
> >> >> >> >>> On Thu, Mar 12, 2009 at 10:25:36AM +0100, Martin Sandve Alnæs wrote:
> >> >> >> >>>> I have serious problems with the idea of letting user-defined
> >> >> >> >>>> functions change nature to discrete functions as a side effect
> >> >> >> >>>> of other operations, effectively hiding the user-defined
> >> >> >> >>>> implementation of eval.
> >> >> >> >>>>
> >> >> >> >>>> (I know this concept wasn't introduced in this discussion, but it's
> >> >> >> >>>> related).
> >> >> >> >>>
> >> >> >> >>> You mean by calling u.vector()? Yes, I agree it's problematic.
> >> >> >> >>>
> >> >> >> >>> But I don't know how to handle it otherwise. Consider the following
> >> >> >> >>> example:
> >> >> >> >>>
> >> >> >> >>>  u = Function(V)
> >> >> >> >>>  solve(A, u.vector(), b)
> >> >> >> >>>
> >> >> >> >>> First u is a user-defined function since it doesn't have a
> >> >> >> >>> vector. Then it becomes discrete when we ask for the vector.
> >> >> >> >>>
> >> >> >> >>> The problem I think is that there is no way for us to check whether a
> >> >> >> >>> user has overloaded eval() without trying to call it.
> >> >> >> >>
> >> >> >> >> If we didn't require that user-defined functions had a function space
> >> >> >> >> for various operations, this wouldn't be as large a problem.
> >> >> >> >> Then you could do
> >> >> >> >>
> >> >> >> >> class MyFunction(Function)
> >> >> >> >> {
> >> >> >> >>   MyFunction(V): Function(V) {}
> >> >> >> >>   MyFunction(): Function() {}
> >> >> >> >>   void eval(...) { ... }
> >> >> >> >> }
> >> >> >> >>
> >> >> >> >> MyFunction f;
> >> >> >> >> f.vector(); // error, no function space!
> >> >> >> >>
> >> >> >> >> MyFunction f(V);
> >> >> >> >> f.vector(); // ok, should interpolate and make f discrete
> >> >> >> >>
> >> >> >> >> This is already possible, but function spaces get attached
> >> >> >> >> unneccesarily:
> >> >> >> >>
> >> >> >> >>
> >> >> >> >> 1) It's not necessary to attach function spaces to functions for
> >> >> >> >> assembly, since a user-defined function is evaluated through
> >> >> >> >> ufc::finite_element::evaluate_dofs for each cell anyway.
> >> >> >> >> This would remove the need for the side-effect in
> >> >> >> >>
> >> >> >> >> MyFunction f;
> >> >> >> >> MyForm a;
> >> >> >> >> a.f = f; // attaches function space to f
> >> >> >
> >> >> > Isn't the function space used to deduct the value rank during assemble now? If
> >> >> > we go for this approach we probably need to readd the rank function and make
> >> >> > it virtual?
> >> >>
> >> >> I don't see how the current approach with rank deducted from
> >> >> some function space is any safer than deducting it from the
> >> >> finite element of the form. Which is to say: this is unsafe, the
> >> >> user should define the rank when defining eval.
> >> >
> >> > We just removed the rank() and dim() arguments from the Function
> >> > interface as they are not needed.
> >>
> >> Are not needed because a Function always has a FunctionSpace
> >> or for some other reason? If we don't require a FunctionSpace
> >> I think they should be re-added. I expect to be able to use a Function
> >> for other things than assembly, e.g., if I get some arbitrary Function
> >> to some code:
> >>
> >> void foo(Function & f)
> >> {
> >>    ... compute something with f
> >> }
> >>
> >> I can't know the dimensions of arguments to f.
> >> This means I can't write generic code to do something with Functions.
> >
> > One of the simplifying assumptions we made a while back when
> > redesigning the Function classes was that a Function always has a
> > FunctionSpace. This still seems to be a good assumption and perhaps we
> > should enforce it even stronger than we do now.
> >
> > It would also make it easy to check for errors during assembly. The
> > form knows the correct FunctionSpace and the assembler may compare
> > that FunctionSpace against the FunctionSpace of the Function.
> 
> I'm happy with "a Function always has a FunctionSpace" concept if it is
> strictly enforced at construction time. The current situation is not so,
> and since I'm allowed to create a user function without a function space
> I expect to be able to use it where it makes sense.
> 
> C++ is often much easier to work with if objects are completely initialized at
> construction time, then there are a lot of mistakes you're just not
> allowed to make.

I agree completely. This is the approach we tried to follow in the
recent redesign of Function. That's why one first needs to create a
mesh, then a function space, then a form. I'm not sure why we
decided not to enforce the requirement of having a function space for
coefficients.

It seems it wouldn't be anymore difficult than what we have now. For
example, the Poisson demo would be

  UnitSquare mesh(32, 32);
  PoissonFunctionSpace V(mesh);

  PoissonBilinearForm a(V, V);
  PoissonLinearForm L(V);
  Source f(V);
  L.f = f; // without side effect

> >> > I don't think having rank(), dim(), geometric_dimension() as mandatory
> >> > in the Function interface would help. I've seen more examples of
> >> > programs breaking from rank() and dim() being implemented incorrectly
> >> > than I've seen programs breaking from eval() being implemented
> >> > incorrectly. And no error-checking can save a user from typing
> >> >
> >> >  values[100] = x[30]
> >> >
> >> > inside eval(), unless we move away from arrays to some bounds-checking
> >> > type for the arguments to eval.
> >>
> >> Ok, it probably isn't a big problem to not have a check for this in assembly.
> >> The main error this kind of check used to catch was mismatch between
> >> function argument order in forms, which isn't a problem anymore.
> >>
> >>
> >> >> > The geometric dimension is proably deducable from other places during
> >> >> > assemble.
> >> >>
> >> >> Of course, but that doesn't guarantee that the Function and the Form
> >> >> uses the same dimensions. Error checking is _the_ reason for
> >> >> making the user provide this information explicitly.
> >> >>
> >> >> The point of having explicit rank in a Function is not that the expected
> >> >> rank isn't available during assembly, but because the actual rank of the
> >> >> Function should be _checked against_ the rank expected for the form argument.
> >> >>
> >> >> In general, I find it very frustrating that DOLFIN has so little
> >> >> error checking.
> >> >
> >> > No one objects to error-checking. Error-checking is good and we should
> >> > do more of it. Feel free to add as much as you like. Make sure to
> >> > include clear and sensible error messages.
> >>
> >> I'm sorry, I don't have any time to donate for this, my time is
> >> better spent with UFL anyway.
> >>
> >> I'm trying to adapt a habit of always documenting my assumptions
> >> of input arguments to a function with checks and assertions.
> >> It's way harder to add checks to other peoples code or code I
> >> wrote myself long ago.
> >>
> >> I believe the only way to do better at this point is if everybody
> >> starts writing checks whenever new code is added or old
> >> code is modified. It also makes reading code easier!
> >
> > Yes, that works fine if you know what you want to do from the start,
> > but I find that the design is often evolving and I wouldn't know which
> > unit tests to write until the code is finished.
> 
> Tests are hard, I know... But here I was thinking about simple
> things like matching dimensions and vector sizes etc.

Yes, and I agree it's important. The tests are missing simply because we've
either been lazy or been pressed on time (and most often the latter).

-- 
Anders

Attachment: signature.asc
Description: Digital signature