dolfin team mailing list archive

[Johan Hake <hake@xxxxxxxxx>]

On Thursday 18 December 2008 00:08:52 Bartosz Sawicki wrote:
> On 17/12/08 05:11 AM, Johan Hake wrote:
> > On Tuesday 16 December 2008 23:30:17 Bartosz Sawicki wrote:
> >> On 16/12/08 03:04 PM, Garth N. Wells wrote:
> >>> Bartosz Sawicki wrote:
> >>>> Before the function revolution there was method cell() inside the
> >>>> Function. I used it to determine index of current cell during assembly
> >>>> procedures, in eval() method. I realized that it disappeared now.
> >>>> How similar functionality can be achieved today?
> >>>
> >>> You'll need to use the eval variant
> >>>
> >>>    void eval(double* values, const Data& data) const;
> >>>
> >>> and you can then access the cell via
> >>>
> >>>   const Cell& cell = data.cell();
> >>
> >> Thanks. That's clear now.
> >> Remember that corresponding section of the manual needs to be updated.
> >
> > That is probably true, together with a lot of other stuff... :)
> >
> >> Do you plan to extend python interface to follow this changes?
> >
> > The pure python interface only support values and x as arguments to the
> > callback function eval(). It should probably be possible to expose the
> > data structure to python too, but this wont be a priority right now.
> >
> > You can always define your own C++ function and compile it in python.
> >
> >    cppcode = '''
> >     class MyFunc : public Function
> >     {
> >     public:
> >
> >       MyFunc(FunctionSpace& V) : Function(V) {}
> >
> >       void eval(double* values, const Data& data)
> >       {
> >          // write your C++ code here!
> >       }
> >
> >     };'''
> >
> >    my_func = Function(V,cppcode)
> >
> > my_func can then be used in form formulation, and then automaticly used
> > in the assembly too. Unfortunatly you cannot use this function directly
> > from the python prompt, i.e.
> >
> >    my_func.eval(v,x)
> >
> > wont work, as the data structure is created duing assemble.
> >
> > Btw, with the new Function interface we have made it much easier to use
> > compiled c++ function from python. Have a look in the docstring of
> > Function. I strongly encourage you to use compiled function instead of
> > pure python functions as the performance boost is significant during
> > assemble.
>
> Thank you for explanation. I've just tried tried to apply this, but I
> failed.
> Please have a look into this simple example:
>
> u = Function("poisson.xml")
> mesh = u.function_space().mesh()
> Vv = VectorFunctionSpace(mesh, "Discontinuous Lagrange", 0)
> cppcode = '''
> class MyFunc : public Function
> {
> public:
>    MyFunc(FunctionSpace& V) : Function(V) {}
>    void eval(double* values, const Data& data){
>      values[0] = 1.0;
>    }
> };
> '''
> v = TestFunction(Vv)
> E = TrialFunction(Vv)
> mat = Function(Vv, cppcode)
> a = dot(v,E)*dx
> L = -mat[0]*dot(v,grad(u))*dx
> pde = LinearPDE(a, L)
> E = pde.solve()
>
> And error message:
>
> Ordering mesh (finished).
> Creating linear PDE with 0 boundary condition(s).
> Solving linear PDE.
>    Assembling matrix over cells (finished).
>    Assembling vector over cells (finished).
> Traceback (most recent call last):
>    File "grad.py", line 41, in <module>
>      E = pde.solve()
>    File
> "/home/sawickib/dolfin/local/lib/python2.5/site-packages/dolfin/pde.py",
> line 54, in solve
>      cpp.LinearPDE.solve(self, u)
>    File
> "/home/sawickib/dolfin/local/lib/python2.5/site-packages/dolfin/cpp.py",
> line 10908, in solve
>      return _cpp.LinearPDE_solve(*args)
> RuntimeError: *** Error: Missing eval() for user-defined function
> (must be overloaded).
>
>
> Error message is a little bit strange for me. Why it doesn't see eval()
> function?

You need to make it const. This is my fault. I never tested that the functions 
in the docstring example actually could be assembled. I will update the 
docstring accordingly. 

Thanks for testing this out for us BArtek!

Johan