Re: [Libmesh-users] Given a physical (or reference, outside of domain) point, which element does it belong to?

[Roy S. <roy...@ic...>]

Oh, and now that I've noticed: please keep Cc:ing the mailing list;
the more helpful discussion that gets archived where future users'
search engines can find it, the better.

Thanks,
---
Roy

On Thu, 16 Nov 2017, Roy Stogner wrote:

>
> On Wed, 15 Nov 2017, Zack Vitoh wrote:
>
>> Pure virtual functions (like those found in the PointLocatorBase and 
>> related classes) were completely new to me earlier today, but I believe I
>> understand the proper syntax, at least, so if it's of use to anyone else, 
>> here is an example of one way to use the PointLocatorBase class (to find
>> the element 'elem_ploc' containing (0,-0.5,0)) 
>> 
>> UniquePtr<PointLocatorBase> 
>> my_point_locator(PointLocatorBase::build(TREE_ELEMENTS, mesh));
>
> This should work, but it's not the most efficient way to go:
> although PointLocatorTree::operator() is O(log N),
> PointLocatorTree::build() is O(N), so you only get a fast amortized
> lookup if you can reuse the same point locator over and over again.
>
> Try
>
> UniquePtr<PointLocatorBase> my_point_locator = mesh.sub_point_locator();
>
> That will create a sub-locator which reuses the same main locator
> instead of building a new one each time.
>
>> Real mpl_tol = 2.0 * diam;
>> my_point_locator->set_close_to_point_tol (mpl_tol);
>> my_point_locator->enable_out_of_mesh_mode();
>
> I assume diam is an element diameter?  Then you're trying to find
> points as far as two diameters away from any current element?  I'm
> afraid that's not guaranteed to work - if you have quads or other
> non-affine elements in your mesh, you can have mapping functions which
> are invertible on the elements (so the mesh is perfectly valid) but
> which become singular far away from the elements (so the
> transformations we do when checking whether an element contains a
> point become invalid).  Beware.
>
> Also, with a huge tolerance, you are going to have multiple elements
> which "contain" a point, and the point locator may not return an
> element which *actually* contains the point, even if one exists, if it
> finds a merely close by element first.
>
>> const Elem* elem_ploc = my_point_locator->operator()( Point(0.,-0.5,0.) );
>
> For operator(), terser syntax is:
>
> const Elem* elem_ploc = (*my_point_locator)( Point(0.,-0.5,0.) );
>
> Knowing the full ugly syntax is still useful, unfortunately, for
> debugging with gdb...
> ---
> Roy