[PPL-devel] [GIT] ppl/ppl(floating_point): Do not use log2 to compute the base 2 logarithm of

[Fabio Bossi]

Module: ppl/ppl
Branch: floating_point
Commit: f6167277a1078a7891e62097d4e034830001f373
URL:    http://www.cs.unipr.it/git/gitweb.cgi?p=ppl/ppl.git;a=commit;h=f6167277a1078a7891e62097d4e034830001f373

Author: Fabio Bossi <bossi at cs.unipr.it>
Date:   Tue Aug  3 11:09:28 2010 +0200

Do not use log2 to compute the base 2 logarithm of
an integer.

---

 src/Float.defs.hh            |    9 +++++++++
 src/Float.inlines.hh         |   33 +++++++++++++++++++++++++++++++++
 src/Float.templates.hh       |    3 +--
 src/Linear_Form.templates.hh |    3 +--
 4 files changed, 44 insertions(+), 4 deletions(-)

diff --git a/src/Float.defs.hh b/src/Float.defs.hh
index 4e66c05..b437184 100644
--- a/src/Float.defs.hh
+++ b/src/Float.defs.hh
@@ -369,6 +369,15 @@ public:
 };
 #endif
 
+// FIXME: is this the right place for this function?
+#ifdef PPL_DOXYGEN_INCLUDE_IMPLEMENTATION_DETAILS
+//! Returns the position of the most significative bit in \p a.
+#endif // defined(PPL_DOXYGEN_INCLUDE_IMPLEMENTATION_DETAILS)
+unsigned int ld2(unsigned long long a);
+
+/* FIXME: some of the following  documentation should probably be
+   under PPL_DOXYGEN_INCLUDE_IMPLEMENTATION_DETAILS */
+
 /*! \brief
   Returns <CODE>true</CODE> if and only if there is some floating point
   number that is representable by \p f2 but not by \p f1.
diff --git a/src/Float.inlines.hh b/src/Float.inlines.hh
index 219797b..2b9c3f2 100644
--- a/src/Float.inlines.hh
+++ b/src/Float.inlines.hh
@@ -451,6 +451,39 @@ is_less_precise_than(Floating_Point_Format f1, Floating_Point_Format f2) {
   return f1 < f2;
 }
 
+#if defined(__GNUC__)
+inline unsigned int ld2(unsigned long long a) {
+ return __builtin_clzll(a) ^ (sizeof(a)*8 - 1);
+}
+#else
+unsigned int ld2(unsigned long long v) {
+ unsigned r = 0;
+ if (v >= 0x100000000ULL) {
+   v >>= 32;
+   r += 32;
+ }
+ if (v >= 0x10000) {
+   v >>= 16;
+   r += 16;
+ }
+ if (v >= 0x100) {
+   v >>= 8;
+   r += 8;
+ }
+ if (v >= 0x10) {
+   v >>= 4;
+   r += 4;
+ }
+ if (v >= 4) {
+   v >>= 2;
+   r += 2;
+ }
+ if (v >= 2)
+   r++;
+ return r;
+}
+#endif
+
 template <typename FP_Interval_Type>
 inline void
 affine_form_image(std::map<dimension_type,
diff --git a/src/Float.templates.hh b/src/Float.templates.hh
index f062fde..287fb4f 100644
--- a/src/Float.templates.hh
+++ b/src/Float.templates.hh
@@ -113,8 +113,7 @@ const FP_Interval_Type& compute_absolute_error(
 
   // We assume that f_base is a power of 2.
   analyzer_format omega;
-  int power = static_cast<int>(log2(f_base)) *
-              (1 - f_exponent_bias - f_mantissa_bits);
+  int power = ld2(f_base) * (1 - f_exponent_bias - f_mantissa_bits);
   omega = std::max(static_cast<analyzer_format>(ldexpl(1.0, power)),
                    std::numeric_limits<analyzer_format>::denorm_min());
 
diff --git a/src/Linear_Form.templates.hh b/src/Linear_Form.templates.hh
index 5bb64c0..596e095 100644
--- a/src/Linear_Form.templates.hh
+++ b/src/Linear_Form.templates.hh
@@ -419,8 +419,7 @@ Linear_Form<C>::relative_error(
 
   C error_propagator;
   // We assume that f_base is a power of 2.
-  int power = static_cast<int>(log2(f_base)) *
-              (-f_mantissa_bits);
+  int power = ld2(f_base) * (-f_mantissa_bits);
   analyzer_format lb = -static_cast<analyzer_format>(ldexpl(1.0, power));
 
   error_propagator.build(i_constraint(GREATER_OR_EQUAL, lb),