7
$\begingroup$

Lets say I've decomposed a continuous function $y(x)$ over some domain $L_x$ using known finite element method with local basis $Q_i(x)$. Suppose $L$ is divided into $M$ "elements". If I want to know the function $y(x)$ at a point $x=p$ (where $p$ is not at a nodal value) then I will need to first know which element $p$ is in. I can then interpolate between the element nodes and find the value $y(p)$. I can do this for all $p$ in $L$ and get back my continuous function.

To do this, one could search through every element to find the where $p$ lies. This is fine in simple cases, but what about more complex cases when there are 1000's of elements and we are working in 3D? Is there an efficient algorithm that I'm missing that could compute this?

$\endgroup$
4
$\begingroup$

Don't do it.

On general meshes, finding an arbitrary point $p$ requires $O(N)$ operations where $N$ is the number of cells. You can get away with significantly less, say $O(\log N)$, if you have a mesh that is constructed hierachically from a coarse mesh by mesh refinement, but unless you use a completely structured grid, finding an arbitrary point inside a mesh is simply a very expensive operation.

Secondly, once you have the cell $K$ in which $p$ lies, you need to find the location $\hat p$ on the reference cell $\hat K$ that corresponds to $p$. For this, you have to invert the mapping $\phi_K: \hat K \mapsto K$, and this (in general) requires a nonlinear iteration -- which is also expensive.

As a consequence, most algorithms try to avoid this at any cost. For example, when you do integration of functions, you loop over all cells, and on each cell over the quadrature points (which are defined on the reference cell, and so only need to be mapped forward with $\phi_K$). This is a far faster approach that looping over a bunch of quadrature points defined in real space (as opposed to reference space relative to each cell).

Ultimately, the point is that you need to re-think the overall approach to doing things so that you always start with a loop over all cells, if at all possible.

$\endgroup$
  • 3
    $\begingroup$ That's great for integration, but what do you suggest when you really do want to know the value of a function a specific point or small set of points? $\endgroup$ – Bill Barth Dec 18 '15 at 17:53
  • $\begingroup$ @BillBarth -- in many cases, that just means that you should thoroughly think about why you need to evaluate the solution there. There are of course legitimate reasons for this, but in most cases, one can recast the algorithm that requires that in ways that yield predictable evaluation points. $\endgroup$ – Wolfgang Bangerth Dec 22 '15 at 2:39
2
$\begingroup$

You could use a domain decomposition approach; divide the domain into coarse subdomains recursively and search through each one (like a binary tree search in 1D).

If you know more about the point $p$, you can probably make it more efficient. For example, if the point $p$ does not change, lookup information can be precomputed and stored during the setup stage. If it changes in time but depends on properties of the solution, it can often be determined by transforming to the reference element and doing the search in parallel.

$\endgroup$
2
$\begingroup$

Similar to comments from Jesse Chan, you could consider breaking up the domain into rectangular (prisms in 3-D) shapes. Then, one option could be to pre-compute a spatial hash table with the finite elements by using a 2-D (or 3-D) coordinate to generate the key. You could then create the table by hashing the vertices of each element.

To efficiently find what element to do computation in, you could hash the coordinates of your arbitrary point $p$ and then loop through the elements in the resulting hash table bucket until you find the element you need to do local computation in. If you build the hash function well, you could end up with very quick look ups.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.