SCOREC / core

parallel finite element unstructured meshes
Other
179 stars 63 forks source link

Expensive gmi_is_point_in_region called during ma::adapt #444

Closed bobpaw closed 3 days ago

bobpaw commented 1 month ago

I'm reporting an issue in ma::adapt where an unnecessarily expensive call (depending on geometric complexity) is made to gmi_is_point_in_region.

Call stack (from gdb):

#0  gmi_is_point_in_region (m=0x8df910, e=0xb, point=0x7ffffffed030)
    at gmi/gmi.c:143
#1  0x00007ffff7cc366d in apf::Mesh::isParamPointInsideModel (this=this@entry=0x7f9d80,
    g=g@entry=0x2a, param=..., x=...)
    at apf/apfMesh.cc:271
#2  0x00007ffff7e77589 in ma::interpolateParametricCoordinatesOnRegularFace (p=..., b=...,
    a=..., t=0.5, g=0x2a, m=0x7f9d80)
    at ma/maSnap.cc:500
#3  ma::interpolateParametricCoordinatesOnFace (p=..., b=..., a=..., t=0.5, g=0x2a,
    m=0x7f9d80) at ma/maSnap.cc:570
#4  ma::interpolateParametricCoordinates (m=0x7f9d80, g=0x2a, t=0.5, a=..., b=..., p=...)
    at ma/maSnap.cc:586
#5  0x00007ffff7e7781b in ma::transferParametricBetween (p=..., t=0.5, v=0x7ffffffed0c0,
    g=0x2a, m=0x7f9d80) at ma/maSnap.cc:601
#6  ma::transferParametricOnEdgeSplit (m=m@entry=0x7f9d80, e=e@entry=0xc9ea2, t=t@entry=0.5,
    p=...) at ma/maSnap.cc:615
#7  0x00007ffff7e5fecb in ma::makeSplitVert (r=r@entry=0x6f6d550, edge=edge@entry=0xc9ea2)
    at ma/maRefine.cc:144
#8  0x00007ffff7e62d8a in ma::splitEdge (r=0x6f6d550, edge=0xc9ea2, v=0x7ffffffed240)
    at ma/maTemplates.cc:81
#9  0x00007ffff7e601af in ma::splitElement (r=r@entry=0x6f6d550, e=e@entry=0xc9ea2)
    at ma/maRefine.cc:244
#10 0x00007ffff7e60643 in ma::splitElements (r=r@entry=0x6f6d550)
    at ma/maRefine.cc:334
#11 0x00007ffff7e60dbd in ma::refine (a=a@entry=0x816270)
    at ma/maRefine.cc:437
#12 0x00007ffff7e5bf77 in ma::adaptVerbose (in=in@entry=0x8e0220, verbose=verbose@entry=false)
    at ma/ma.cc:79
#13 0x000000000040450e in main (argc=<optimized out>, argv=<optimized out>)
    at capstone_clis/capAdapt.cc:972

I'm not particularly well-versed in the adapt code, but this should be avoidable during adaptation.

cwsmith commented 1 month ago

I think there are scenarios where it makes sense to run the check. For example, if an mesh edge classified on the model boundary is being split then the point where it is split could be outside the domain of the geometric model. A crude drawing is attached. splitEdge

joshia5 commented 1 month ago

for reference, can you add details on the adaptation inputs you're using? or a the link to the main function if available

bobpaw commented 1 month ago

for reference, can you add details on the adaptation inputs you're using? or a the link to the main function if available

This are the inputs we're using.

ma::Input* in;
in = ma::makeAdvanced(ma::configure(apfMesh, mdsScaleField, mdsFrameField));
in->shouldSnap = true;
in->shouldTransferParametric = true;
in->shouldFixShape = true;
in->shouldForceAdaptation = true;
if (apfMesh->getDimension() == 2)
  in->goodQuality = 0.04; // this is mean-ratio squared
 else // 3D meshes
   in->goodQuality = 0.027; // this is the mean-ratio cubed
in->maximumIterations = 10;

We are calling it from a file similar to capstone_clis/capStoneAnisoAdapt.cc

cwsmith commented 1 month ago

For reference, this commit https://github.com/SCOREC/core/commit/4b514947ae542ea3a0bb436deac041d32cc0e343 added the check.

bobpaw commented 3 days ago

I believe this is closed by #446?

cwsmith commented 3 days ago

Yeah. It should be all set.