Subject: [FFC-dev] vertexeval function

From: "Garth N. Wells" <[email-masked]>
Date: Mon, 29 May 2006 19:10:23 +0200

For evaluating functions at vertxes, FFC writes a function

void vertexeval(real values[], unsigned int vertex, const real x[],
const Mesh& mesh) const

Would there be a problem in changing "const real x[]" to "const
GenericVector x"? This would make it easier to represent a function
using different vector data structures.

Garth

---

From: Anders Logg <[email-masked]>
Date: Mon, 29 May 2006 19:32:42 +0200

> using different vector data structures.

The problem would be that vertexeval() would then need to call
GenericMatrix::operator() (uint i) which is slow for PETScSparseMatrix.

But perhaps calling PETScSparseMatrix::array() is just as slow?

vertexeval() is just a temporary fix for Lagrange elements. Somehow,
we need to compute the values at all vertices so we can write the
solution to file. What we should really do would be to compute the
projection onto piecewise linears and then just pick the values.

Maybe we could always generate code for the projection onto piecewise
linears when a form is compiled in FFC? If U is in V' and V is the
space of piecewise linears, then we need to compute the projection
W of U onto V:

    int v W dx = int v U dx for all v in V

On matrix form, we have

    MW = b

where b represents the right-hand side int U v dx.

The mass matrix M is already computed by default in DOLFIN so we don't
need to generate it. We only need to generate the linear form v*U*dx.

Then we could remove vertexeval.

/Anders

---

From: "Garth N. Wells" <[email-masked]>
Date: Mon, 29 May 2006 19:57:35 +0200

> But perhaps calling PETScSparseMatrix::array() is just as slow?

When manipulating a lot of term, PETScSparseMatrix::array() is much
faster. I think that this has to with the calls

  VecAssemblyBegin(x);
  VecAssemblyEnd(x);

after accessing each entry.

> projection onto piecewise linears and then just pick the values.

OK. A simple temporary fix would be if vextexeval returned an array of
integers containing the positions in the big vector. Then the values
could be pulled as a block out of the vector. 

> need to generate it. We only need to generate the linear form v*U*dx.

. . . and solve the linear system.

Garth