diff --git a/Packages/Interval_arithmetic/changes.txt b/Packages/Interval_arithmetic/changes.txt index ed2cc4b90f8..69711620b20 100644 --- a/Packages/Interval_arithmetic/changes.txt +++ b/Packages/Interval_arithmetic/changes.txt @@ -1,5 +1,8 @@ Changes done to the Interval Arithmetic package. +Version 4.84 on 17 April 2001 +- New filtered predicate : coplanar_side_of_oriented_circleC3. + Version 4.83 on 10 April 2001 - FPU.h is changed to reflect the realtity on Alpha/Linux 2.2.18 with libc 2.1 to support FPE handling by the Linux kernel. diff --git a/Packages/Interval_arithmetic/include/CGAL/Arithmetic_filter/predicates/kernel_ftC3.h b/Packages/Interval_arithmetic/include/CGAL/Arithmetic_filter/predicates/kernel_ftC3.h index 6c67dcf9a16..1716468970f 100644 --- a/Packages/Interval_arithmetic/include/CGAL/Arithmetic_filter/predicates/kernel_ftC3.h +++ b/Packages/Interval_arithmetic/include/CGAL/Arithmetic_filter/predicates/kernel_ftC3.h @@ -1058,6 +1058,358 @@ orientationC3( #endif // CGAL_IA_NEW_FILTERS +#ifndef CGAL_CFG_MATCHING_BUG_2 +template < class CGAL_IA_CT, class CGAL_IA_ET, bool CGAL_IA_PROTECTED, + class CGAL_IA_CACHE > +#else +static +#endif +/* CGAL_KERNEL_LARGE_INLINE */ +Oriented_side +coplanar_side_of_oriented_circleC3( + const Filtered_exact &px, + const Filtered_exact &py, + const Filtered_exact &pz, + const Filtered_exact &qx, + const Filtered_exact &qy, + const Filtered_exact &qz, + const Filtered_exact &rx, + const Filtered_exact &ry, + const Filtered_exact &rz, + const Filtered_exact &tx, + const Filtered_exact &ty, + const Filtered_exact &tz) +{ + try + { + Protect_FPU_rounding Protection; + return coplanar_side_of_oriented_circleC3( + px.interval(), + py.interval(), + pz.interval(), + qx.interval(), + qy.interval(), + qz.interval(), + rx.interval(), + ry.interval(), + rz.interval(), + tx.interval(), + ty.interval(), + tz.interval()); + } + catch (Interval_nt_advanced::unsafe_comparison) + { + Protect_FPU_rounding Protection(CGAL_FE_TONEAREST); + return coplanar_side_of_oriented_circleC3( + px.exact(), + py.exact(), + pz.exact(), + qx.exact(), + qy.exact(), + qz.exact(), + rx.exact(), + ry.exact(), + rz.exact(), + tx.exact(), + ty.exact(), + tz.exact()); + } +} + +#ifdef CGAL_IA_NEW_FILTERS + +struct Static_Filtered_coplanar_side_of_oriented_circleC3_12 +{ + static double _bound; + static double _epsilon_0; + static unsigned number_of_failures; // ? + static unsigned number_of_updates; + + static Oriented_side update_epsilon( + const Static_filter_error &px, + const Static_filter_error &py, + const Static_filter_error &pz, + const Static_filter_error &qx, + const Static_filter_error &qy, + const Static_filter_error &qz, + const Static_filter_error &rx, + const Static_filter_error &ry, + const Static_filter_error &rz, + const Static_filter_error &tx, + const Static_filter_error &ty, + const Static_filter_error &tz, + double & epsilon_0) + { + typedef Static_filter_error FT; + + + + + + FT ptx = px - tx; + FT pty = py - ty; + FT ptz = pz - tz; + FT pt2 = CGAL_NTS square(ptx) + CGAL_NTS square(pty) + CGAL_NTS square(ptz); + FT qtx = qx - tx; + FT qty = qy - ty; + FT qtz = qz - tz; + FT qt2 = CGAL_NTS square(qtx) + CGAL_NTS square(qty) + CGAL_NTS square(qtz); + FT rtx = rx - tx; + FT rty = ry - ty; + FT rtz = rz - tz; + FT rt2 = CGAL_NTS square(rtx) + CGAL_NTS square(rty) + CGAL_NTS square(rtz); + FT pqx = qx - px; + FT pqy = qy - py; + FT pqz = qz - pz; + FT prx = rx - px; + FT pry = ry - py; + FT prz = rz - pz; + FT vx = pqy*prz - pqz*pry; + FT vy = pqz*prx - pqx*prz; + FT vz = pqx*pry - pqy*prx; + FT v2 = CGAL_NTS square(vx) + CGAL_NTS square(vy) + CGAL_NTS square(vz); + return Oriented_side(Static_Filtered_sign_of_determinant4x4_16::update_epsilon(ptx,pty,ptz,pt2, + rtx,rty,rtz,rt2, + qtx,qty,qtz,qt2, + vx,vy,vz,v2, + epsilon_0)); + } + + // Call this function from the outside to update the context. + static void new_bound (const double b) // , const double error = 0) + { + _bound = b; + number_of_updates++; + // recompute the epsilons: "just" call it over Static_filter_error. + // That's the tricky part that might not work for everything. + (void) update_epsilon(b,b,b,b,b,b,b,b,b,b,b,b,_epsilon_0); + // TODO: We should verify that all epsilons have really been updated. + } + + static Oriented_side epsilon_variant( + const Restricted_double &px, + const Restricted_double &py, + const Restricted_double &pz, + const Restricted_double &qx, + const Restricted_double &qy, + const Restricted_double &qz, + const Restricted_double &rx, + const Restricted_double &ry, + const Restricted_double &rz, + const Restricted_double &tx, + const Restricted_double &ty, + const Restricted_double &tz, + const double & epsilon_0) + { + typedef Restricted_double FT; + + + + + + FT ptx = px - tx; + FT pty = py - ty; + FT ptz = pz - tz; + FT pt2 = CGAL_NTS square(ptx) + CGAL_NTS square(pty) + CGAL_NTS square(ptz); + FT qtx = qx - tx; + FT qty = qy - ty; + FT qtz = qz - tz; + FT qt2 = CGAL_NTS square(qtx) + CGAL_NTS square(qty) + CGAL_NTS square(qtz); + FT rtx = rx - tx; + FT rty = ry - ty; + FT rtz = rz - tz; + FT rt2 = CGAL_NTS square(rtx) + CGAL_NTS square(rty) + CGAL_NTS square(rtz); + FT pqx = qx - px; + FT pqy = qy - py; + FT pqz = qz - pz; + FT prx = rx - px; + FT pry = ry - py; + FT prz = rz - pz; + FT vx = pqy*prz - pqz*pry; + FT vy = pqz*prx - pqx*prz; + FT vz = pqx*pry - pqy*prx; + FT v2 = CGAL_NTS square(vx) + CGAL_NTS square(vy) + CGAL_NTS square(vz); + return Oriented_side(Static_Filtered_sign_of_determinant4x4_16::epsilon_variant(ptx,pty,ptz,pt2, + rtx,rty,rtz,rt2, + qtx,qty,qtz,qt2, + vx,vy,vz,v2, + epsilon_0)); + } +}; + +#ifndef CGAL_CFG_MATCHING_BUG_2 +template < class CGAL_IA_CT, class CGAL_IA_ET, class CGAL_IA_CACHE > +#else +static +#endif +/* CGAL_KERNEL_LARGE_INLINE */ +Oriented_side +coplanar_side_of_oriented_circleC3( + const Filtered_exact &px, + const Filtered_exact &py, + const Filtered_exact &pz, + const Filtered_exact &qx, + const Filtered_exact &qy, + const Filtered_exact &qz, + const Filtered_exact &rx, + const Filtered_exact &ry, + const Filtered_exact &rz, + const Filtered_exact &tx, + const Filtered_exact &ty, + const Filtered_exact &tz) +{ +// bool re_adjusted = false; + const double SAF_bound = Static_Filtered_coplanar_side_of_oriented_circleC3_12::_bound; + + // Check the bounds. All arguments must be <= SAF_bound. + if ( + fabs(px.to_double()) > SAF_bound || + fabs(py.to_double()) > SAF_bound || + fabs(pz.to_double()) > SAF_bound || + fabs(qx.to_double()) > SAF_bound || + fabs(qy.to_double()) > SAF_bound || + fabs(qz.to_double()) > SAF_bound || + fabs(rx.to_double()) > SAF_bound || + fabs(ry.to_double()) > SAF_bound || + fabs(rz.to_double()) > SAF_bound || + fabs(tx.to_double()) > SAF_bound || + fabs(ty.to_double()) > SAF_bound || + fabs(tz.to_double()) > SAF_bound) + { +// re_adjust: + // Compute the new bound. + double NEW_bound = 0.0; + NEW_bound = max(NEW_bound, fabs(px.to_double())); + NEW_bound = max(NEW_bound, fabs(py.to_double())); + NEW_bound = max(NEW_bound, fabs(pz.to_double())); + NEW_bound = max(NEW_bound, fabs(qx.to_double())); + NEW_bound = max(NEW_bound, fabs(qy.to_double())); + NEW_bound = max(NEW_bound, fabs(qz.to_double())); + NEW_bound = max(NEW_bound, fabs(rx.to_double())); + NEW_bound = max(NEW_bound, fabs(ry.to_double())); + NEW_bound = max(NEW_bound, fabs(rz.to_double())); + NEW_bound = max(NEW_bound, fabs(tx.to_double())); + NEW_bound = max(NEW_bound, fabs(ty.to_double())); + NEW_bound = max(NEW_bound, fabs(tz.to_double())); + // Re-adjust the context. + Static_Filtered_coplanar_side_of_oriented_circleC3_12::new_bound(NEW_bound); + } + + try + { + return Static_Filtered_coplanar_side_of_oriented_circleC3_12::epsilon_variant( + px.dbl(), + py.dbl(), + pz.dbl(), + qx.dbl(), + qy.dbl(), + qz.dbl(), + rx.dbl(), + ry.dbl(), + rz.dbl(), + tx.dbl(), + ty.dbl(), + tz.dbl(), + Static_Filtered_coplanar_side_of_oriented_circleC3_12::_epsilon_0); + } + catch (...) + { + // if (!re_adjusted) { // It failed, we re-adjust once. + // re_adjusted = true; + // goto re_adjust; + // } + Static_Filtered_coplanar_side_of_oriented_circleC3_12::number_of_failures++; + return coplanar_side_of_oriented_circleC3( + px.exact(), + py.exact(), + pz.exact(), + qx.exact(), + qy.exact(), + qz.exact(), + rx.exact(), + ry.exact(), + rz.exact(), + tx.exact(), + ty.exact(), + tz.exact()); + } +} + +#ifndef CGAL_CFG_MATCHING_BUG_2 +template < class CGAL_IA_CT, class CGAL_IA_ET, class CGAL_IA_CACHE > +#else +static +#endif +/* CGAL_KERNEL_LARGE_INLINE */ +Oriented_side +coplanar_side_of_oriented_circleC3( + const Filtered_exact &px, + const Filtered_exact &py, + const Filtered_exact &pz, + const Filtered_exact &qx, + const Filtered_exact &qy, + const Filtered_exact &qz, + const Filtered_exact &rx, + const Filtered_exact &ry, + const Filtered_exact &rz, + const Filtered_exact &tx, + const Filtered_exact &ty, + const Filtered_exact &tz) +{ + CGAL_assertion_code( + const double SAF_bound = Static_Filtered_coplanar_side_of_oriented_circleC3_12::_bound; ) + CGAL_assertion(!( + fabs(px.to_double()) > SAF_bound || + fabs(py.to_double()) > SAF_bound || + fabs(pz.to_double()) > SAF_bound || + fabs(qx.to_double()) > SAF_bound || + fabs(qy.to_double()) > SAF_bound || + fabs(qz.to_double()) > SAF_bound || + fabs(rx.to_double()) > SAF_bound || + fabs(ry.to_double()) > SAF_bound || + fabs(rz.to_double()) > SAF_bound || + fabs(tx.to_double()) > SAF_bound || + fabs(ty.to_double()) > SAF_bound || + fabs(tz.to_double()) > SAF_bound)); + + try + { + return Static_Filtered_coplanar_side_of_oriented_circleC3_12::epsilon_variant( + px.dbl(), + py.dbl(), + pz.dbl(), + qx.dbl(), + qy.dbl(), + qz.dbl(), + rx.dbl(), + ry.dbl(), + rz.dbl(), + tx.dbl(), + ty.dbl(), + tz.dbl(), + Static_Filtered_coplanar_side_of_oriented_circleC3_12::_epsilon_0); + } + catch (...) + { + Static_Filtered_coplanar_side_of_oriented_circleC3_12::number_of_failures++; + return coplanar_side_of_oriented_circleC3( + px.exact(), + py.exact(), + pz.exact(), + qx.exact(), + qy.exact(), + qz.exact(), + rx.exact(), + ry.exact(), + rz.exact(), + tx.exact(), + ty.exact(), + tz.exact()); + } +} + +#endif // CGAL_IA_NEW_FILTERS + #ifndef CGAL_CFG_MATCHING_BUG_2 template < class CGAL_IA_CT, class CGAL_IA_ET, bool CGAL_IA_PROTECTED, class CGAL_IA_CACHE > diff --git a/Packages/Interval_arithmetic/include/CGAL/Arithmetic_filter/static_infos/predicates/kernel_ftC3.h b/Packages/Interval_arithmetic/include/CGAL/Arithmetic_filter/static_infos/predicates/kernel_ftC3.h index bddc674f1e3..04de3569e56 100644 --- a/Packages/Interval_arithmetic/include/CGAL/Arithmetic_filter/static_infos/predicates/kernel_ftC3.h +++ b/Packages/Interval_arithmetic/include/CGAL/Arithmetic_filter/static_infos/predicates/kernel_ftC3.h @@ -58,6 +58,13 @@ unsigned Static_Filtered_orientationC3_12::number_of_updates = 0; unsigned Static_Filtered_orientationC3_12::number_of_failures = 0; +double Static_Filtered_coplanar_side_of_oriented_circleC3_12::_epsilon_0; +double Static_Filtered_coplanar_side_of_oriented_circleC3_12::_bound = -1.0; + +unsigned Static_Filtered_coplanar_side_of_oriented_circleC3_12::number_of_updates = 0; + +unsigned Static_Filtered_coplanar_side_of_oriented_circleC3_12::number_of_failures = 0; + double Static_Filtered_equal_directionC3_6::_epsilon_0; double Static_Filtered_equal_directionC3_6::_epsilon_1; double Static_Filtered_equal_directionC3_6::_epsilon_2;