On Oct 16, 2009, at 6:04 PM, Anders Logg wrote:
On Fri, Oct 16, 2009 at 03:34:06PM -0400, mspieg wrote:
Hi All,
I'm trying to update a working 0.9.3 code to 0.9.4, and I need
to understand
the new interfaces for Function::eval a bit better.
I had things working fine with the old interface
Function::eval (double *values, const double *x, const ufc::cell
&ufc_cell,
uint cell_index)
and am trying to understand how to efficiently use
Function::eval (double *values, const double *x, const Cell
&dolfin_cell, const
ufc::cell &ufc_cell)
(My codes need to do a large number of evals at arbitrary points
and cells
efficiently)
among other things, is there an easy way to update the index of a
cell, without
constructing a new Cell every time? (or shouldn't I be worrying
about this).
Constructing a cell should not cost that much. It just stores the
index, the dimension and a reference to the mesh.
If you are using iterators over a mesh, then you already have the cell
by dereferencing the pointer. The typical usage would be
UFCCell ufc_cell(mesh);
for (CellIterator cell(mesh); !mesh.end(); ++cell)
{
ufc_cell.update(*cell);
f.eval(values, x, *cell, ufc_cell);
}
But perhaps you don't use iterators?
Thanks Anders,
I haven't been using iterators because I've been overloading eval
(double* values, Data& data) and getting the cells passed through
data during assembly
the issue is that I'm not necessarily evaluating at the position and
cell in data, but at another point and cell derived from the current
location and cell (it's basically a particle tracking problem).
Ideally, given data I would set another Data structure (data_prime)
with the new point, cell and ufc_cell, then call eval as needed. As
long as the cost to construct a cell is small, it will probably be
easy enough to do with the current code (the problem is I need to do
this at every quadrature point, and I found previously that
construction of UFCCell at every point was exceptionally expensive).
Anyway, I'll play around a bit more and if I have specific issues or
design questions, I'll let you know.
cheers
marc
Anyway, I need to play with things a bit more, but looking at the
source code
for eval, I was also wondering
if lines 320 and 321 are required as line 322 doesn't seem to
require cell.
typo? bug? or general ignorance on my part (probably the latter).
Yes, lines 320-321 should be removed. Thanks for noticing. I've fixed
this now.
--
Anders
cheers
marc
void Function::eval(double* values, const Data& data) const
00311 {
00312 assert(values);
00313 assert(_function_space);
00314
00315 // Check if UFC cell if available and cell matches
00316 if (data._dolfin_cell && _function_space->has_cell
(*data._dolfin_cell))
00317 {
00318 // Efficient evaluation on given cell
00319 assert(data._ufc_cell);
00320 const uint cell_index = data._ufc_cell->entity_indices
[data._ufc_cell->topological_dimension][0];
00321 Cell cell(_function_space->mesh(), cell_index);
00322 eval(values, data.x, *data._dolfin_cell, *data._ufc_cell);
00323 }
00324 else
00325 {
00326 // Redirect to point-based evaluation
00327 eval(values, data.x);
00328 }
00329 }
----------------------------------------------------
Marc Spiegelman
Lamont-Doherty Earth Observatory
Dept. of Applied Physics/Applied Math
Columbia University
http://www.ldeo.columbia.edu/~mspieg
tel: 845 704 2323 (SkypeIn)
----------------------------------------------------
_______________________________________________
DOLFIN-dev mailing list
DOLFIN-dev@fenics.org
http://www.fenics.org/mailman/listinfo/dolfin-dev
_______________________________________________
DOLFIN-dev mailing list
DOLFIN-dev@fenics.org
http://www.fenics.org/mailman/listinfo/dolfin-dev
----------------------------------------------------
Marc Spiegelman
Lamont-Doherty Earth Observatory
Dept. of Applied Physics/Applied Math
Columbia University
http://www.ldeo.columbia.edu/~mspieg
tel: 845 704 2323 (SkypeIn)
----------------------------------------------------
_______________________________________________
DOLFIN-dev mailing list
DOLFIN-dev@fenics.org
http://www.fenics.org/mailman/listinfo/dolfin-dev
