← Back to team overview

dolfin team mailing list archive

Re: Evaluating the FEM solution at an arbitrary point

 

On Thu, Feb 14, 2008 at 03:44:38PM +0100, Dag Lindbo wrote:
> Anders Logg wrote:
> > On Thu, Feb 14, 2008 at 12:00:38PM +0100, Dag Lindbo wrote:
> >>> Then I suggest first finding out which cells those points lie in, then
> >>> then for each cell with a point get the expansion coefficients within
> >>> that cell, then multiply those coefficients with the values of the
> >>> basis functions at the points.
> >>>
> >>> The basis functions are available from the ufc::finite_element.
> >>>
> >> The searching part has been sorted out. How do I get the 
> >> ufc::finite_element from the function and/or the form? E.g.
> >>
> >> Function f(mesh, 0.0);
> >> StokesTHBilinearForm a;
> >> StokesTHLinearForm L(f);
> >> LinearPDE pde(a, L, mesh, bcs);
> >>
> >> pde.set("PDE linear solver", "direct");
> >> pde.solve(u, p);
> >>
> >> ufc::form frm = a.form(); // form() is virtual
> >> ufc::finite_element fe = frm.create_finite_element();
> > 
> > You should be able to do this:
> > 
> >   // Create element (the '0' means create it for the space used for the
> >   // first argument in the form, '1' for second etc)
> >   ufc::finite_element* element = a.form().create_finite_element(0);
> 
> Thanks! In my case element->space_dimension() return 15 (dofs) which 
> makes sense, since I have a Taylor-Hood element (but will problably 
> switch later). However, I only want to pick out the velocity components. 
> How should a mixed element be handled?

This is what the functions create_sub_element and create_sub_dof_map
are for. If you have an element (or dof_map), you can pick out the
element for any sub systems. You can have arbitrary depths of sub
system nesting. For example in your case, you have

  w = (u, p)

so you can ask for the finite element for w and then ask that element
for the elements for u and p. Then again, u is a system of size d
(space dimension) so then you can again ask for the sub systems for
the individual components.

Look in the UFC manual for details on the create_sub_foo functions.

-- 
Anders


> I have(u is a Function):
> 
>    Cell cell(mesh,cell_nr);
>    UFCCell ufc_cell(cell);
> 
>    ufc::finite_element* fe = a->form().create_finite_element(0);
>    ufc::dof_map* ufc_dof_map = a->form().create_dof_map(0);
> 
>    DofMap* dof_map = new DofMap(*ufc_dof_map, mesh);
> 
>    unsigned int* dofs = new unsigned int[dof_map->local_dimension()];
>    real* basis_values = new real[fe->space_dimension()];
>    real* dof_values = new real[fe->space_dimension()];
> 
>    Vector vec = u.vector();
> 
>    dof_map->tabulate_dofs(dofs,ufc_cell);
>    vec.get(dof_values, dof_map->local_dimension(), dofs);
> 
> > 
> >   // Use it for a while
> >   element->tabulate_basis(...);
> >   ...
> > 
> >   // Then delete it
> >   delete element;
> > 
> _______________________________________________
> DOLFIN-dev mailing list
> DOLFIN-dev@xxxxxxxxxx
> http://www.fenics.org/mailman/listinfo/dolfin-dev


Follow ups

References