@@ -966,28 +966,6 @@ static bool triangle_area_is_zero(const ccd_vec3_t& a, const ccd_vec3_t& b,
966
966
return false ;
967
967
}
968
968
969
- /* *
970
- * Determines if the point P is on the line segment AB.
971
- * If A, B, P are coincident, report true.
972
- */
973
- static bool is_point_on_line_segment (const ccd_vec3_t & p, const ccd_vec3_t & a,
974
- const ccd_vec3_t & b) {
975
- if (are_coincident (a, b)) {
976
- return are_coincident (a, p);
977
- }
978
- // A and B are not coincident, if the triangle ABP has non-zero area, then P
979
- // is not on the line adjoining AB, and hence not on the line segment AB.
980
- if (!triangle_area_is_zero (a, b, p)) {
981
- return false ;
982
- }
983
- // P is on the line adjoinging AB. If P is on the line segment AB, then
984
- // PA.dot(PB) <= 0.
985
- ccd_vec3_t PA, PB;
986
- ccdVec3Sub2 (&PA, &p, &a);
987
- ccdVec3Sub2 (&PB, &p, &b);
988
- return ccdVec3Dot (&PA, &PB) <= 0 ;
989
- }
990
-
991
969
/* *
992
970
* Computes the normal vector of a triangular face on a polytope, and the normal
993
971
* vector points outward from the polytope. Notice we assume that the origin
@@ -1691,6 +1669,8 @@ static int __ccdGJK(const void *obj1, const void *obj2,
1691
1669
*/
1692
1670
static void validateNearestFeatureOfPolytopeBeingEdge (ccd_pt_t * polytope) {
1693
1671
assert (polytope->nearest_type == CCD_PT_EDGE);
1672
+
1673
+ constexpr ccd_real_t kEps = constants<ccd_real_t >::eps ();
1694
1674
// Only verify the feature if the nearest feature is an edge.
1695
1675
1696
1676
const ccd_pt_edge_t * const nearest_edge =
@@ -1709,27 +1689,29 @@ static void validateNearestFeatureOfPolytopeBeingEdge(ccd_pt_t* polytope) {
1709
1689
// n̂ ⋅ (o - vₑ) ≤ 0 or, with simplification, -n̂ ⋅ vₑ ≤ 0 (since n̂ ⋅ o = 0).
1710
1690
origin_to_face_distance[i] =
1711
1691
-ccdVec3Dot (&face_normals[i], &nearest_edge->vertex [0 ]->v .v );
1712
- if (origin_to_face_distance[i] > 0 ) {
1692
+ // If the origin lies *on* the edge, then it also lies on the two adjacent
1693
+ // faces. Rather than failing on strictly *positive* signed distance, we
1694
+ // introduce an epsilon threshold. This usage of epsilon is to account for a
1695
+ // discrepancy in the signed distance computation. How GJK (and partially
1696
+ // EPA) compute the signed distance from origin to face may *not* be exactly
1697
+ // the same as done in this test (i.e. for a given set of vertices, the
1698
+ // plane equation can be defined various ways. Those ways are
1699
+ // *mathematically* equivalent but numerically will differ due to rounding).
1700
+ // We account for those differences by allowing a very small positive signed
1701
+ // distance to be considered zero. We assume that the GJK/EPA code
1702
+ // ultimately clasifies inside/outside around *zero* and *not* epsilon.
1703
+ if (origin_to_face_distance[i] > kEps ) {
1713
1704
FCL_THROW_FAILED_AT_THIS_CONFIGURATION (
1714
1705
" The origin is outside of the polytope. This should already have "
1715
1706
" been identified as separating." );
1716
1707
}
1717
1708
}
1718
- // We compute the projection of the origin onto the plane as
1719
- // -face_normals[i] * origin_to_face_distance[i]
1720
- // If the projection to both faces are on the edge, the the edge is the
1721
- // closest feature.
1722
- bool is_edge_closest_feature = true ;
1723
- for (int i = 0 ; i < 2 ; ++i) {
1724
- ccd_vec3_t origin_projection_to_plane = face_normals[i];
1725
- ccdVec3Scale (&(origin_projection_to_plane), -origin_to_face_distance[i]);
1726
- if (!is_point_on_line_segment (origin_projection_to_plane,
1727
- nearest_edge->vertex [0 ]->v .v ,
1728
- nearest_edge->vertex [1 ]->v .v )) {
1729
- is_edge_closest_feature = false ;
1730
- break ;
1731
- }
1732
- }
1709
+
1710
+ // We know the reported squared distance to the edge. If that distance is
1711
+ // functionally zero, then the edge must *truly* be the nearest feature.
1712
+ // If it isn't, then it must be one of the adjacent faces.
1713
+ const bool is_edge_closest_feature = nearest_edge->dist < kEps * kEps ;
1714
+
1733
1715
if (!is_edge_closest_feature) {
1734
1716
// We assume that libccd is not crazily wrong. Although the closest
1735
1717
// feature is not the edge, it is near that edge. Hence we select the
@@ -1779,7 +1761,7 @@ static int __ccdEPA(const void *obj1, const void *obj2,
1779
1761
return -2 ;
1780
1762
}
1781
1763
1782
- while (1 ){
1764
+ while (1 ) {
1783
1765
// get triangle nearest to origin
1784
1766
*nearest = ccdPtNearest (polytope);
1785
1767
if (polytope->nearest_type == CCD_PT_EDGE) {
@@ -1791,14 +1773,13 @@ static int __ccdEPA(const void *obj1, const void *obj2,
1791
1773
*nearest = ccdPtNearest (polytope);
1792
1774
}
1793
1775
1794
- // get next support point
1795
- if (nextSupport (polytope, obj1, obj2, ccd, *nearest, &supp) != 0 ) {
1796
- break ;
1797
- }
1776
+ // get next support point
1777
+ if (nextSupport (polytope, obj1, obj2, ccd, *nearest, &supp) != 0 ) {
1778
+ break ;
1779
+ }
1798
1780
1799
- // expand nearest triangle using new point - supp
1800
- if (expandPolytope (polytope, *nearest, &supp) != 0 )
1801
- return -2 ;
1781
+ // expand nearest triangle using new point - supp
1782
+ if (expandPolytope (polytope, *nearest, &supp) != 0 ) return -2 ;
1802
1783
}
1803
1784
1804
1785
return 0 ;
0 commit comments