awslabs / palace

3D finite element solver for computational electromagnetics
https://awslabs.github.io/palace/dev
Apache License 2.0
224 stars 50 forks source link

Fix bounding ball calculation for coaxial lumped ports #246

Closed sebastiangrimberg closed 1 month ago

sebastiangrimberg commented 1 month ago

Previously, for coaxial lumped ports with particularly coarse meshes, the following assertion would fail in geodata.cpp's BoundingBallFromPointCloud:

MFEM_VERIFY(std::abs(PerpendicularDistance({delta}, perp) - ball.radius) <=
                rel_tol * ball.radius,
            "Furthest point perpendicular must be on the exterior of the ball: "
                << PerpendicularDistance({delta}, perp) << " vs. " << ball.radius
                << "!");

This PR improves the bounding ball algorithm to avoid this error (and is a simpler algorithm). We now just compute the ball origin as the average of all points and radial direction as the point furthest from the centroid. The rest of the algorithm is unchanged.

Here is an example second-order mesh which was causing problems prior to this PR:

Screenshot 2024-05-06 at 4 14 13 PM

sebastiangrimberg commented 1 month ago

Converted to draft while revising the bounding box/bounding ball algorithms.

sebastiangrimberg commented 1 month ago

Commit c22bf04 corrects the previous commit after noticing that the oriented bounding box calculation from mesh::GetBoundingBox only provides a correct result when the convex hull of the points forms a rectangle or rectangular prism. For spheres/circles, the box is inscribed in the circle/sphere so the use in CoaxialElementData needs adjusting.

This could be an opportunity to implement a true oriented bounding box calculation, if it might be useful elsewhere. One option is from Baraquet and Har-Peled, where the core ingredient would be a 2D minimal bounding rectangle implementation. See, for example, the implementation here which is Boost Licensed.

Another idea is to first compute the convex hull of the points, followed by an SVD-based approach for computing the approximate oriented bounding box.

Alternatively, for now, the current implementation works fine for boundaries which are rectangles and circles as we expect in current usage.

sebastiangrimberg commented 1 month ago

Following up on this comment, the latest PR implements Welzl's algorithm to compute the minimum bounding sphere which is used for coaxial ports now. The comment for when the simplified OBB calculation given by mesh::ComputeBoundingBox has been clarified.