Firstly: Very good idea.  We have our own mesh to mesh projection code here 
(mostly for interpolating high order solutions on fine low order grids for 
visualization)... and it does take a bit to run for 3D solutions... and the 
speedup for periodic boundary conditions and some of our point source code will 
definitely help the runtime of some of our applications.

As for the actual problem... I don't think that there is any calculation that 
you can do with just the nodes that will _guarantee_ that all of the element is 
contained in the bounding box.  I think the only foolproof way to do it is to 
use the gradients of the shape functions to actually find the minima / maxima 
(this will be a small newton solve I believe).  Performing that operation for 
each side in 3D will probably be more expensive than just doing inverse_map() 
in the first place....

Now... that is the theory.  However, I think that in practice our elements are 
MUCH more well behaved than the pathological cases.  For instance, it's pretty 
hard for me to imagine an element being so distorted that it lies outside of a 
bounding sphere made using the corner nodes.  In fact, I kind of wonder if such 
an element would even have a positive jacobian....

Does this mean we should have a "secure = true" option for contains_point() 
just like we do in other places?  I would vote for when "secure = false" we 
build a bounding box and add a fudge factor based on the relative distance 
between any two nodes.  I think that in practice this will save a lot of time 
and be correct 99.99% of the time....

Derek

On Jun 15, 2011, at 11:58 AM, Roy Stogner wrote:

> 
> I just caught myself nearly introducing a potential bug into
> libMesh... but the code that would trigger this corner-case bug is a
> *really* tempting change, so I'm wondering if anyone has ideas on how
> to fix it.
> 
> 
> Elem::contains_point() calls FE::inverse_map().  But FE::inverse_map()
> is expensive, and most of the time in this context it's unnecessary -
> the tested point is well outside the element and we don't care exactly
> where outside.
> 
> So my thought was to precede the inverse_map() with a bounding box
> check: if the point tested has a higher x value than any of the
> element's nodes (+ TOLERANCE), for example, then contains_point() can
> return false and there's no need to bother with more expensive tests.
> 
> Trying this on a mesh-to-mesh projection code of mine (which I'll add
> to src/apps if people want it, by the way - our grid2grid app was
> hard-coded for older I/O formats...) gave me an *8-fold* speedup in
> opt mode.  That was just on a 2D case; the savings in 3D are likely
> greater.  The speedup might be noticeable for anyone doing periodic
> boundary conditions or other PointLocator based operations too.
> 
> 
> Here's the problem: how do we get such a test to avoid false negatives
> on elements with curved (i.e. potentially bulging out of the bounding
> box) sides?  I've got the svn head using a bounding box test only on
> first-order elements, which should work (first-order hexes can have
> curved sides, but the curves never extend outside of the bounding
> box), and that's better than nothing, but I'd like a test (even an
> alternate less discriminating test) applicable to higher-order
> geometric elements.
> ---
> Roy
> 
> ------------------------------------------------------------------------------
> EditLive Enterprise is the world's most technically advanced content
> authoring tool. Experience the power of Track Changes, Inline Image
> Editing and ensure content is compliant with Accessibility Checking.
> http://p.sf.net/sfu/ephox-dev2dev
> _______________________________________________
> Libmesh-devel mailing list
> [email protected]
> https://lists.sourceforge.net/lists/listinfo/libmesh-devel


------------------------------------------------------------------------------
EditLive Enterprise is the world's most technically advanced content
authoring tool. Experience the power of Track Changes, Inline Image
Editing and ensure content is compliant with Accessibility Checking.
http://p.sf.net/sfu/ephox-dev2dev
_______________________________________________
Libmesh-devel mailing list
[email protected]
https://lists.sourceforge.net/lists/listinfo/libmesh-devel

Reply via email to