## Re: [Libmesh-devel] Patch

 Re: [Libmesh-devel] Patch From: Roy Stogner - 2006-08-07 20:18:39 ```On Mon, 7 Aug 2006, Tim Kr=F6ger wrote: > On Fri, 4 Aug 2006, Roy Stogner wrote: > >> The inverse map iteration definitely won't converge for some elements >> if the physical point is outside the element, but that shouldn't cause >> a problem for Elem::constains_point() anymore; if its call to >> inverse_map() doesn't converge, the result is a distant point that >> will make sure on_reference_element() is false. > > No: If the nonlinear solver doesn't converge, it will call error() and = thus=20 > crash. At least, this call to error() should be removed -- perhaps by=20 > giving the solver class an optional argument. That's exactly why the inverse_map function has the "bool secure" argument. If secure is false (as it is when contains_point() calls inverse_map() ), then Newton isn't expected to converge and inverse_map just exits after 10 iterations. Of course, inverse_map() is *supposed* to exit by returning an out-of-bounds point, not just whatever its last iterate was. Derek, didn't we fix this? Or did you send me a patch I forgot to commit to CVS? > The other problem that I see when leaving the linear search in is that = when=20 > I will later implement the possibiliy to evaluate the MeshFunction outs= ide=20 > the domain covered by the grid, the linear search will slow this down=20 > considerabely. That is undesirable. > What about adding some heuristic safety margin around the bounding box = in=20 > the case of quadratic mappings? Well, not really satisfying ... As long as the margin is a proveable upper bound, I'll be happy with it. Trying to figure out what would be an upper bound on a quadratic Hex makes my head spin, though. > Perhaps better: Enable the linear search only if quadratic mappings are= =20 > involved. That sounds good too, but deciding whether or not quadratic mappings are involved is easier said than done. I'd go so far as to say that second-order geometric elements with linear mappings is one of the most common usage patterns in libMesh, since we need the extra nodes in the second-order elements to store degree of freedom numbers for any p>1 finite element. > Any suggestions? No good ones. All I can think to do is warn about the problem. When you add an interface for out-of-mesh MeshFunction extensions, let that interface turn off (well, replace) the linear search, and add an assert() or #ifdef DEBUG block in that code path which verifies that it doesn't see any quadratic mappings. --- Roy```