dolfin team mailing list archive

["Martin Sandve Alnæs" <martinal@xxxxxxxxx>]

2008/11/3 Anders Logg <logg@xxxxxxxxx>:
> On Mon, Nov 03, 2008 at 03:55:31PM +0000, Garth N. Wells wrote:
>>
>>
>> Anders Logg wrote:
>> > On Mon, Nov 03, 2008 at 02:38:36PM +0000, Garth N. Wells wrote:
>> >>
>> >> Anders Logg wrote:
>> >>> On Mon, Nov 03, 2008 at 11:22:26AM +0000, Garth N. Wells wrote:
>> >>>> Anders Logg wrote:
>> >>>>> On Sun, Nov 02, 2008 at 06:29:25PM +0000, Garth N. Wells wrote:
>> >>>>>> Anders Logg wrote:
>> >>>>>>> On Sun, Nov 02, 2008 at 05:52:21PM +0000, Garth N. Wells wrote:
>> >>>>>>>> Do we want to insist that Dirichlet bc functions that do not appear
>> >>>>>>>> inside a form are constructed with a FunctionSpace? DirichletBC is
>> >>>>>>>> supplied with a FunctionSpace, so if the bc Function does not have a
>> >>>>>>>> FunctionSpace, we could attach one automatically.
>> >>>>>>>>
>> >>>>>>>> Garth
>> >>>>>>>> _______________________________________________
>> >>>>>>>> DOLFIN-dev mailing list
>> >>>>>>>> DOLFIN-dev@xxxxxxxxxx
>> >>>>>>>> http://www.fenics.org/mailman/listinfo/dolfin-dev
>> >>>>>>> I think this is already handled. Look in the Poisson demo. It uses a
>> >>>>>>> Constant to set the BC and it does not have a FunctionSpace attached
>> >>>>>>> to it. The DirichletBC class now uses its own FunctionSpace rather
>> >>>>>>> than the one that the Function has (if any). There is a check (in
>> >>>>>>> DirichletBC::check()) that checks that the FunctionSpace for the
>> >>>>>>> Function is the same as the one in the DirichletBC.
>> >>>>>>>
>> >>>>>> It works for Constant, but not for Functions. I was getting an error
>> >>>>>> when Function::interpolate is called. Function::interpolate leads to
>> >>>>>> eval being called, in which case there is a test for the FunctionSpace
>> >>>>>> which fails. Constant provides its own eval and therefore doesn't have a
>> >>>>>> problem.
>> >>>>>>
>> >>>>>> For now, I've added a test in DirichletBC for the FunctionSpace. What we
>> >>>>>> can add is an attach function if there is no FunctionSpace associated.
>> >>>>>>
>> >>>>>> Garth
>> >>>>> In which demo does this show up? Is there a simple way I can comment
>> >>>>> something out to reproduce the error so I understand what goes wrong?
>> >>>>>
>> >>>> Look at /demo/nls/nonlinearpoisson/cpp.
>> >>>>
>> >>>> If you change
>> >>>>
>> >>>>    DirichletBoundaryCondition g(V, t);
>> >>>>
>> >>>> to
>> >>>>
>> >>>>    DirichletBoundaryCondition g(t);
>> >>>>
>> >>>> it will break down.
>> >>>>
>> >>>> Garth
>> >>> ok I see the problem now.
>> >>>
>> >>> The problem is a user may choose to either overload a scalar eval
>> >>> function or a tensor eval function and we need to decide which one
>> >>> after the callback from ufc::function::evaluate(). If the
>> >>> FunctionSpace is not known, we can't decide which one to pick.
>> >>>
>> >>> If we insist that one should be able to pass a Function without a
>> >>> FunctionSpace to a DirichletBC, then we must remove the scalar eval
>> >>> function.
>> >>>
>> >> Fine with me. I think that it makes things simpler because the eval
>> >> interface remains the same for all user-defined functions.
>> >>
>> >> Garth
>> >
>> > ok. It will also look the same as in Python.
>> >
>>
>> We discussed recently passing an object to eval() which contains some
>> data. It would be useful the object also carried information on the rank
>> and dimension of the function to allow checks and switching between
>> 1D/2D/3D problems.
>>
>> Garth
>
> Yes, this would be nice, but it won't work as long as we don't require
> that a Function always has a FunctionSpace. We could add some new
> classes to handle error checking and data for user-defined Functions:
>
>  void eval(Values& values, Data& data)
>  {
>    values[0] = sin(data.x[0]);
>  }
>
> The class Values could check that data is not assigned to any illegal
> indices and it could also check that all values have been assigned
> etc, but if the Function does not know its FunctionSpace, this can't
> be done.
>
> Another complication related to this but also to the thread-safety of
> cell() and facet() is that UFC gets in the middle of the call
> sequence:
>
>  assemble()
>  |--> Function::interpolate()
>       |--> FiniteElement::evaluate_dof()
>            |--> ufc::finite_element::evaluate_dof
>                 |--> ufc::function::evaluate()
>                      |--> Function::eval()
>
> Since eval() is called from the generated UFC code, any arguments like
> cell and facet passed from the assembler will be lost on the way.
>
> Should we extend the UFC interface to allow sending a void* to
> evaluate_dof which it will propagate to evaluate()?

That's possible, but not very safe.

A typesafe alternative is that dolfin::Function doesn't inherit ufc::function,
but a ufc::function subclass is created which has a dolfin::Function
which it calls and a dolfin::Data & which it sends in the call.
Then we avoid ufc complications in dolfin::Function, and one such
wrapper function can be created for each thread.
It doesn't really add any more function calls, since it replaces the
existing ufc::function::evaluate in the call stack you wrote above.

dolfin_ufc_function->function = dolfin_function;
for(cell_iterator = ...; ...)
{
  dolfin_ufc_function->data._cell = *cell_iterator;
  ...
  ufc_fe->evaluate_dofs(..., dolfin_ufc_function, ...);
}

--
Martin