[Boost-commit] svn:boost r74501 - in sandbox/multiprecision/boost: . utility

Subject: [Boost-commit] svn:boost r74501 - in sandbox/multiprecision/boost: . utility
From: e_float_at_[hidden]
Date: 2011-09-21 17:00:18

Author: christopher_kormanyos
Date: 2011-09-21 17:00:15 EDT (Wed, 21 Sep 2011)
New Revision: 74501
URL: http://svn.boost.org/trac/boost/changeset/74501

Log:
- Added contents of Boost.Multiprecision.
Added:
   sandbox/multiprecision/boost/   (props changed)
   sandbox/multiprecision/boost/mp_complex.hpp   (contents, props changed)
   sandbox/multiprecision/boost/mp_float.hpp   (contents, props changed)
   sandbox/multiprecision/boost/mp_float_base.hpp   (contents, props changed)
   sandbox/multiprecision/boost/mp_float_efx.hpp   (contents, props changed)
   sandbox/multiprecision/boost/mp_float_functions.hpp   (contents, props changed)
   sandbox/multiprecision/boost/mp_float_gmp.hpp   (contents, props changed)
   sandbox/multiprecision/boost/mp_float_mpfr.hpp   (contents, props changed)
   sandbox/multiprecision/boost/utility/   (props changed)
   sandbox/multiprecision/boost/utility/util_alternating_sum.hpp   (contents, props changed)
   sandbox/multiprecision/boost/utility/util_coefficient_expansion.hpp   (contents, props changed)
   sandbox/multiprecision/boost/utility/util_digit_scale.hpp   (contents, props changed)
   sandbox/multiprecision/boost/utility/util_find_root_base.hpp   (contents, props changed)
   sandbox/multiprecision/boost/utility/util_find_root_bisect.hpp   (contents, props changed)
   sandbox/multiprecision/boost/utility/util_find_root_newton_raphson.hpp   (contents, props changed)
   sandbox/multiprecision/boost/utility/util_function_base.hpp   (contents, props changed)
   sandbox/multiprecision/boost/utility/util_function_derivative.hpp   (contents, props changed)
   sandbox/multiprecision/boost/utility/util_function_operation.hpp   (contents, props changed)
   sandbox/multiprecision/boost/utility/util_interpolate.hpp   (contents, props changed)
   sandbox/multiprecision/boost/utility/util_numeric_cast.hpp   (contents, props changed)
   sandbox/multiprecision/boost/utility/util_point.hpp   (contents, props changed)
   sandbox/multiprecision/boost/utility/util_power_j_pow_x.hpp   (contents, props changed)
   sandbox/multiprecision/boost/utility/util_power_x_pow_n.hpp   (contents, props changed)
   sandbox/multiprecision/boost/utility/util_ranged_function_operation.hpp   (contents, props changed)
   sandbox/multiprecision/boost/utility/util_timer.hpp   (contents, props changed)
   sandbox/multiprecision/boost/utility/util_trapezoid.hpp   (contents, props changed)

Added: sandbox/multiprecision/boost/mp_complex.hpp
==============================================================================
--- (empty file)
+++ sandbox/multiprecision/boost/mp_complex.hpp	2011-09-21 17:00:15 EDT (Wed, 21 Sep 2011)
@@ -0,0 +1,497 @@
+
+//          Copyright Christopher Kormanyos 2002 - 2011.
+// Distributed under the Boost Software License, Version 1.0.
+//    (See accompanying file LICENSE_1_0.txt or copy at
+//          http://www.boost.org/LICENSE_1_0.txt)
+
+// This work is based on an earlier work:
+// "Algorithm 910: A Portable C++ Multiple-Precision System for Special-Function Calculations",
+// in ACM TOMS, {VOL 37, ISSUE 4, (February 2011)} (C) ACM, 2011. http://doi.acm.org/10.1145/1916461.1916469
+
+#ifndef _MP_COMPLEX_HPP_
+  #define _MP_COMPLEX_HPP_
+
+  #include <complex>
+
+  #include <boost/multiprecision/mp_float.hpp>
+
+  // A separate complex class for mp_float has been created. Even though
+  // a generic template class std::complex<T> exists, the C++ specification
+  // ISO7IEC 14882:2003 paragraph 26.2/2 indicates that: "The effect of
+  // instantiating the template complex<T> for any type other than float,
+  // double or long double is unspecified". The strict interpretation thereof
+  // disallows both using the template class std::complex<T> with mp_float
+  // as well as creating a template specialization for mp_float.
+  // Therefore the separate class mp_complex is needed.
+
+  namespace boost
+  {
+    namespace multiprecision
+    {
+      class mp_complex;
+
+      inline mp_complex& operator+=(mp_complex&, const char);
+      inline mp_complex& operator+=(mp_complex&, const signed char);
+      inline mp_complex& operator+=(mp_complex&, const signed short);
+      inline mp_complex& operator+=(mp_complex&, const signed int);
+      inline mp_complex& operator+=(mp_complex&, const signed long);
+      inline mp_complex& operator+=(mp_complex&, const signed long long);
+      inline mp_complex& operator+=(mp_complex&, const unsigned char);
+      inline mp_complex& operator+=(mp_complex&, const wchar_t);
+      inline mp_complex& operator+=(mp_complex&, const unsigned short);
+      inline mp_complex& operator+=(mp_complex&, const unsigned int);
+      inline mp_complex& operator+=(mp_complex&, const unsigned long);
+      inline mp_complex& operator+=(mp_complex&, const unsigned long long);
+      inline mp_complex& operator+=(mp_complex&, const float);
+      inline mp_complex& operator+=(mp_complex&, const double);
+      inline mp_complex& operator+=(mp_complex&, const long double);
+
+      inline mp_complex& operator-=(mp_complex&, const char);
+      inline mp_complex& operator-=(mp_complex&, const signed char);
+      inline mp_complex& operator-=(mp_complex&, const signed short);
+      inline mp_complex& operator-=(mp_complex&, const signed int);
+      inline mp_complex& operator-=(mp_complex&, const signed long);
+      inline mp_complex& operator-=(mp_complex&, const signed long long);
+      inline mp_complex& operator-=(mp_complex&, const unsigned char);
+      inline mp_complex& operator-=(mp_complex&, const wchar_t);
+      inline mp_complex& operator-=(mp_complex&, const unsigned short);
+      inline mp_complex& operator-=(mp_complex&, const unsigned int);
+      inline mp_complex& operator-=(mp_complex&, const unsigned long);
+      inline mp_complex& operator-=(mp_complex&, const unsigned long long);
+      inline mp_complex& operator-=(mp_complex&, const float);
+      inline mp_complex& operator-=(mp_complex&, const double);
+      inline mp_complex& operator-=(mp_complex&, const long double);
+
+      inline mp_complex& operator*=(mp_complex&, const char);
+      inline mp_complex& operator*=(mp_complex&, const signed char);
+      inline mp_complex& operator*=(mp_complex&, const signed short);
+      inline mp_complex& operator*=(mp_complex&, const signed int);
+      inline mp_complex& operator*=(mp_complex&, const signed long);
+      inline mp_complex& operator*=(mp_complex&, const signed long long);
+      inline mp_complex& operator*=(mp_complex&, const unsigned char);
+      inline mp_complex& operator*=(mp_complex&, const wchar_t);
+      inline mp_complex& operator*=(mp_complex&, const unsigned short);
+      inline mp_complex& operator*=(mp_complex&, const unsigned int);
+      inline mp_complex& operator*=(mp_complex&, const unsigned long);
+      inline mp_complex& operator*=(mp_complex&, const unsigned long long);
+      inline mp_complex& operator*=(mp_complex&, const float);
+      inline mp_complex& operator*=(mp_complex&, const double);
+      inline mp_complex& operator*=(mp_complex&, const long double);
+
+      inline mp_complex& operator/=(mp_complex&, const char);
+      inline mp_complex& operator/=(mp_complex&, const signed char);
+      inline mp_complex& operator/=(mp_complex&, const signed short);
+      inline mp_complex& operator/=(mp_complex&, const signed int);
+      inline mp_complex& operator/=(mp_complex&, const signed long);
+      inline mp_complex& operator/=(mp_complex&, const signed long long);
+      inline mp_complex& operator/=(mp_complex&, const unsigned char);
+      inline mp_complex& operator/=(mp_complex&, const wchar_t);
+      inline mp_complex& operator/=(mp_complex&, const unsigned short);
+      inline mp_complex& operator/=(mp_complex&, const unsigned int);
+      inline mp_complex& operator/=(mp_complex&, const unsigned long);
+      inline mp_complex& operator/=(mp_complex&, const unsigned long long);
+      inline mp_complex& operator/=(mp_complex&, const float);
+      inline mp_complex& operator/=(mp_complex&, const double);
+      inline mp_complex& operator/=(mp_complex&, const long double);
+
+      inline mp_complex& operator+=(mp_complex&, const mp_float&);
+      inline mp_complex& operator-=(mp_complex&, const mp_float&);
+      inline mp_complex& operator*=(mp_complex&, const mp_float&);
+      inline mp_complex& operator/=(mp_complex&, const mp_float&);
+    }
+  }
+
+  namespace boost
+  {
+    namespace multiprecision
+    {
+      class mp_complex
+      {
+      private:
+        mp_float Re;
+        mp_float Im;
+
+      public:
+        mp_complex(const char n)               : Re(n), Im(boost::multiprecision::zero()) { }
+        mp_complex(const signed char n)        : Re(n), Im(boost::multiprecision::zero()) { }
+        mp_complex(const unsigned char n)      : Re(n), Im(boost::multiprecision::zero()) { }
+        mp_complex(const wchar_t n)            : Re(n), Im(boost::multiprecision::zero()) { }
+        mp_complex(const signed short n)       : Re(n), Im(boost::multiprecision::zero()) { }
+        mp_complex(const unsigned short n)     : Re(n), Im(boost::multiprecision::zero()) { }
+        mp_complex(const signed int n)         : Re(n), Im(boost::multiprecision::zero()) { }
+        mp_complex(const unsigned int n)       : Re(n), Im(boost::multiprecision::zero()) { }
+        mp_complex(const signed long n)        : Re(n), Im(boost::multiprecision::zero()) { }
+        mp_complex(const unsigned long n)      : Re(n), Im(boost::multiprecision::zero()) { }
+        mp_complex(const signed long long n)   : Re(n), Im(boost::multiprecision::zero()) { }
+        mp_complex(const unsigned long long n) : Re(n), Im(boost::multiprecision::zero()) { }
+        mp_complex(const float f)              : Re(f), Im(boost::multiprecision::zero()) { }
+        mp_complex(const double d)             : Re(d), Im(boost::multiprecision::zero()) { }
+        mp_complex(const long double ld)       : Re(ld), Im(boost::multiprecision::zero()) { }
+
+        explicit mp_complex(const char* const s)    : Re(s), Im(boost::multiprecision::zero()) { }
+        explicit mp_complex(const std::string& str) : Re(str), Im(boost::multiprecision::zero()) { }
+
+        mp_complex() : Re(boost::multiprecision::zero()), Im(boost::multiprecision::zero()) { }
+        mp_complex(const mp_float& re) : Re(re), Im(boost::multiprecision::zero()) { }
+        mp_complex(const mp_float& re, const mp_float& im) : Re(re), Im(im) { }
+        mp_complex(const mp_complex& z) : Re(z.Re), Im(z.Im) { }
+
+        mp_float real(void) const { return Re; }
+        mp_float imag(void) const { return Im; }
+
+        static mp_float real(const mp_complex& z) { return z.Re; }
+        static mp_float imag(const mp_complex& z) { return z.Im; }
+
+        mp_float norm(void) const { return (Re * Re) + (Im * Im); }
+
+        mp_complex& operator=(const mp_complex& v) { Re = v.Re; Im = v.Im; return *this; }
+        mp_complex& operator=(const mp_float& v)    { Re = v; Im = boost::multiprecision::zero(); return *this; }
+
+        mp_complex& operator+=(const mp_complex& v) { Re += v.Re; Im += v.Im; return *this; }
+        mp_complex& operator-=(const mp_complex& v) { Re -= v.Re; Im -= v.Im; return *this; }
+
+        mp_complex& operator*=(const mp_complex& v)
+        {
+          const mp_float tmp_re = (Re * v.Re) - (Im * v.Im);
+          const mp_float tmp_im = (Re * v.Im) + (Im * v.Re);
+          Re = tmp_re;
+          Im = tmp_im;
+          return *this;
+        }
+
+        mp_complex& operator/=(const mp_complex& v)
+        {
+          const mp_float one_over_denom = boost::multiprecision::one() / v.norm();
+          const mp_float tmp_re = ((Re * v.Re) + (Im * v.Im)) * one_over_denom;
+          const mp_float tmp_im = ((Im * v.Re) - (Re * v.Im)) * one_over_denom;
+          Re = tmp_re;
+          Im = tmp_im;
+          return *this;
+        }
+
+        // Operators pre-increment and pre-decrement
+        const mp_complex& operator++(void) { ++Re; return *this; }
+        const mp_complex& operator--(void) { --Re; return *this; }
+
+        bool isnan   (void) const { return Re.isnan()    || Im.isnan(); }
+        bool isinf   (void) const { return Re.isinf()    || Im.isinf(); }
+        bool isfinite(void) const { return Re.isfinite() && Im.isfinite(); }
+        bool isneg   (void) const { return Re.isneg(); }
+        bool ispos   (void) const { return Re.ispos(); }
+        bool isint   (void) const { return Re.isint()  && Im.iszero(); }
+        bool isone   (void) const { return Re.isone()  && Im.iszero(); }
+        bool iszero  (void) const { return Re.iszero() && Im.iszero(); }
+
+        friend inline mp_complex& operator+=(mp_complex&, const char);
+        friend inline mp_complex& operator+=(mp_complex&, const signed char);
+        friend inline mp_complex& operator+=(mp_complex&, const signed short);
+        friend inline mp_complex& operator+=(mp_complex&, const signed int);
+        friend inline mp_complex& operator+=(mp_complex&, const signed long);
+        friend inline mp_complex& operator+=(mp_complex&, const signed long long);
+        friend inline mp_complex& operator+=(mp_complex&, const unsigned char);
+        friend inline mp_complex& operator+=(mp_complex&, const wchar_t);
+        friend inline mp_complex& operator+=(mp_complex&, const unsigned short);
+        friend inline mp_complex& operator+=(mp_complex&, const unsigned int);
+        friend inline mp_complex& operator+=(mp_complex&, const unsigned long);
+        friend inline mp_complex& operator+=(mp_complex&, const unsigned long long);
+        friend inline mp_complex& operator+=(mp_complex&, const float);
+        friend inline mp_complex& operator+=(mp_complex&, const double);
+        friend inline mp_complex& operator+=(mp_complex&, const long double);
+
+        friend inline mp_complex& operator-=(mp_complex&, const char);
+        friend inline mp_complex& operator-=(mp_complex&, const signed char);
+        friend inline mp_complex& operator-=(mp_complex&, const signed short);
+        friend inline mp_complex& operator-=(mp_complex&, const signed int);
+        friend inline mp_complex& operator-=(mp_complex&, const signed long);
+        friend inline mp_complex& operator-=(mp_complex&, const signed long long);
+        friend inline mp_complex& operator-=(mp_complex&, const unsigned char);
+        friend inline mp_complex& operator-=(mp_complex&, const wchar_t);
+        friend inline mp_complex& operator-=(mp_complex&, const unsigned short);
+        friend inline mp_complex& operator-=(mp_complex&, const unsigned int);
+        friend inline mp_complex& operator-=(mp_complex&, const unsigned long);
+        friend inline mp_complex& operator-=(mp_complex&, const unsigned long long);
+        friend inline mp_complex& operator-=(mp_complex&, const float);
+        friend inline mp_complex& operator-=(mp_complex&, const double);
+        friend inline mp_complex& operator-=(mp_complex&, const long double);
+
+        friend inline mp_complex& operator*=(mp_complex&, const char);
+        friend inline mp_complex& operator*=(mp_complex&, const signed char);
+        friend inline mp_complex& operator*=(mp_complex&, const signed short);
+        friend inline mp_complex& operator*=(mp_complex&, const signed int);
+        friend inline mp_complex& operator*=(mp_complex&, const signed long);
+        friend inline mp_complex& operator*=(mp_complex&, const signed long long);
+        friend inline mp_complex& operator*=(mp_complex&, const unsigned char);
+        friend inline mp_complex& operator*=(mp_complex&, const wchar_t);
+        friend inline mp_complex& operator*=(mp_complex&, const unsigned short);
+        friend inline mp_complex& operator*=(mp_complex&, const unsigned int);
+        friend inline mp_complex& operator*=(mp_complex&, const unsigned long);
+        friend inline mp_complex& operator*=(mp_complex&, const unsigned long long);
+        friend inline mp_complex& operator*=(mp_complex&, const float);
+        friend inline mp_complex& operator*=(mp_complex&, const double);
+        friend inline mp_complex& operator*=(mp_complex&, const long double);
+
+        friend inline mp_complex& operator/=(mp_complex&, const char);
+        friend inline mp_complex& operator/=(mp_complex&, const signed char);
+        friend inline mp_complex& operator/=(mp_complex&, const signed short);
+        friend inline mp_complex& operator/=(mp_complex&, const signed int);
+        friend inline mp_complex& operator/=(mp_complex&, const signed long);
+        friend inline mp_complex& operator/=(mp_complex&, const signed long long);
+        friend inline mp_complex& operator/=(mp_complex&, const unsigned char);
+        friend inline mp_complex& operator/=(mp_complex&, const wchar_t);
+        friend inline mp_complex& operator/=(mp_complex&, const unsigned short);
+        friend inline mp_complex& operator/=(mp_complex&, const unsigned int);
+        friend inline mp_complex& operator/=(mp_complex&, const unsigned long);
+        friend inline mp_complex& operator/=(mp_complex&, const unsigned long long);
+        friend inline mp_complex& operator/=(mp_complex&, const float);
+        friend inline mp_complex& operator/=(mp_complex&, const double);
+        friend inline mp_complex& operator/=(mp_complex&, const long double);
+
+        friend inline mp_complex& operator+=(mp_complex&, const mp_float&);
+        friend inline mp_complex& operator-=(mp_complex&, const mp_float&);
+        friend inline mp_complex& operator*=(mp_complex&, const mp_float&);
+        friend inline mp_complex& operator/=(mp_complex&, const mp_float&);
+      };
+    }
+  }
+
+  namespace boost
+  {
+    namespace multiprecision
+    {
+      inline std::ostream& operator<<(std::ostream& os, const mp_complex& z) { return os << '(' << z.real() << ',' << z.imag() << ')'; }
+
+      // Global unary operators of mp_float reference.
+      inline       mp_complex  operator-(const mp_complex& u) { return mp_complex(-u.real(), -u.imag()); }
+      inline       mp_complex& operator+(      mp_complex& u) { return u; }
+      inline const mp_complex& operator+(const mp_complex& u) { return u; }
+
+      // Global operators post-increment and post-decrement
+      inline mp_complex operator++(mp_complex& u, int) { const mp_complex v(u); ++u; return v; }
+      inline mp_complex operator--(mp_complex& u, int) { const mp_complex v(u); --u; return v; }
+
+      // Global comparison operators
+      inline bool operator==(const mp_complex& u, const mp_complex& v) { return (u.real() == v.real()) && (u.imag() == v.imag()); }
+      inline bool operator!=(const mp_complex& u, const mp_complex& v) { return (u.real() != v.real()) || (u.imag() != v.imag()); }
+
+      // Global arithmetic operators with const mp_complex& and const mp_complex&.
+      inline mp_complex operator+(const mp_complex& u, const mp_complex& v) { return mp_complex(u) += v; }
+      inline mp_complex operator-(const mp_complex& u, const mp_complex& v) { return mp_complex(u) -= v; }
+      inline mp_complex operator*(const mp_complex& u, const mp_complex& v) { return mp_complex(u) *= v; }
+      inline mp_complex operator/(const mp_complex& u, const mp_complex& v) { return mp_complex(u) /= v; }
+
+      // Global arithmetic operators with const mp_complex& and const mp_float&.
+      inline mp_complex operator+(const mp_complex& u, const mp_float& v) { return mp_complex(u.real() + v, u.imag()); }
+      inline mp_complex operator-(const mp_complex& u, const mp_float& v) { return mp_complex(u.real() - v, u.imag()); }
+      inline mp_complex operator*(const mp_complex& u, const mp_float& v) { return mp_complex(u.real() * v, u.imag() * v); }
+      inline mp_complex operator/(const mp_complex& u, const mp_float& v) { return mp_complex(u.real() / v, u.imag() / v); }
+
+      // Global arithmetic operators with const mp_float& and const mp_complex&.
+      inline mp_complex operator+(const mp_float& u, const mp_complex& v) { return mp_complex(u + v.real(), v.imag()); }
+      inline mp_complex operator-(const mp_float& u, const mp_complex& v) { return mp_complex(u - v.real(), -v.imag()); }
+      inline mp_complex operator*(const mp_float& u, const mp_complex& v) { return mp_complex(u * v.real(), u * v.imag()); }
+      inline mp_complex operator/(const mp_float& u, const mp_complex& v) { const mp_float v_norm = v.norm(); return mp_complex((u * v.real()) / v_norm, (-u * v.imag()) / v_norm); }
+
+      // Global add/sub/mul/div of const mp_complex& with all built-in types.
+      inline mp_complex operator+(const mp_complex& z, const char n)               { return mp_complex(z.real() + n, z.imag()); }
+      inline mp_complex operator+(const mp_complex& z, const signed char n)        { return mp_complex(z.real() + n, z.imag()); }
+      inline mp_complex operator+(const mp_complex& z, const signed short n)       { return mp_complex(z.real() + n, z.imag()); }
+      inline mp_complex operator+(const mp_complex& z, const signed int n)         { return mp_complex(z.real() + n, z.imag()); }
+      inline mp_complex operator+(const mp_complex& z, const signed long n)        { return mp_complex(z.real() + n, z.imag()); }
+      inline mp_complex operator+(const mp_complex& z, const signed long long n)   { return mp_complex(z.real() + n, z.imag()); }
+      inline mp_complex operator+(const mp_complex& z, const unsigned char n)      { return mp_complex(z.real() + n, z.imag()); }
+      inline mp_complex operator+(const mp_complex& z, const wchar_t n)            { return mp_complex(z.real() + n, z.imag()); }
+      inline mp_complex operator+(const mp_complex& z, const unsigned short n)     { return mp_complex(z.real() + n, z.imag()); }
+      inline mp_complex operator+(const mp_complex& z, const unsigned int n)       { return mp_complex(z.real() + n, z.imag()); }
+      inline mp_complex operator+(const mp_complex& z, const unsigned long n)      { return mp_complex(z.real() + n, z.imag()); }
+      inline mp_complex operator+(const mp_complex& z, const unsigned long long n) { return mp_complex(z.real() + n, z.imag()); }
+      inline mp_complex operator+(const mp_complex& z, const float f)              { return mp_complex(z.real() + f, z.imag()); }
+      inline mp_complex operator+(const mp_complex& z, const double d)             { return mp_complex(z.real() + d, z.imag()); }
+      inline mp_complex operator+(const mp_complex& z, const long double ld)       { return mp_complex(z.real() + ld, z.imag()); }
+
+      inline mp_complex operator-(const mp_complex& z, const char n)               { return mp_complex(z.real() - n, z.imag()); }
+      inline mp_complex operator-(const mp_complex& z, const signed char n)        { return mp_complex(z.real() - n, z.imag()); }
+      inline mp_complex operator-(const mp_complex& z, const signed short n)       { return mp_complex(z.real() - n, z.imag()); }
+      inline mp_complex operator-(const mp_complex& z, const signed int n)         { return mp_complex(z.real() - n, z.imag()); }
+      inline mp_complex operator-(const mp_complex& z, const signed long n)        { return mp_complex(z.real() - n, z.imag()); }
+      inline mp_complex operator-(const mp_complex& z, const signed long long n)   { return mp_complex(z.real() - n, z.imag()); }
+      inline mp_complex operator-(const mp_complex& z, const unsigned char n)      { return mp_complex(z.real() - n, z.imag()); }
+      inline mp_complex operator-(const mp_complex& z, const wchar_t n)            { return mp_complex(z.real() - n, z.imag()); }
+      inline mp_complex operator-(const mp_complex& z, const unsigned short n)     { return mp_complex(z.real() - n, z.imag()); }
+      inline mp_complex operator-(const mp_complex& z, const unsigned int n)       { return mp_complex(z.real() - n, z.imag()); }
+      inline mp_complex operator-(const mp_complex& z, const unsigned long n)      { return mp_complex(z.real() - n, z.imag()); }
+      inline mp_complex operator-(const mp_complex& z, const unsigned long long n) { return mp_complex(z.real() - n, z.imag()); }
+      inline mp_complex operator-(const mp_complex& z, const float f)              { return mp_complex(z.real() - f, z.imag()); }
+      inline mp_complex operator-(const mp_complex& z, const double d)             { return mp_complex(z.real() - d, z.imag()); }
+      inline mp_complex operator-(const mp_complex& z, const long double ld)       { return mp_complex(z.real() - ld, z.imag()); }
+
+      inline mp_complex operator*(const mp_complex& z, const char n)               { return mp_complex(z.real() * n, z.imag() * n); }
+      inline mp_complex operator*(const mp_complex& z, const signed char n)        { return mp_complex(z.real() * n, z.imag() * n); }
+      inline mp_complex operator*(const mp_complex& z, const signed short n)       { return mp_complex(z.real() * n, z.imag() * n); }
+      inline mp_complex operator*(const mp_complex& z, const signed int n)         { return mp_complex(z.real() * n, z.imag() * n); }
+      inline mp_complex operator*(const mp_complex& z, const signed long n)        { return mp_complex(z.real() * n, z.imag() * n); }
+      inline mp_complex operator*(const mp_complex& z, const signed long long n)   { return mp_complex(z.real() * n, z.imag() * n); }
+      inline mp_complex operator*(const mp_complex& z, const unsigned char n)      { return mp_complex(z.real() * n, z.imag() * n); }
+      inline mp_complex operator*(const mp_complex& z, const wchar_t n)            { return mp_complex(z.real() * n, z.imag() * n); }
+      inline mp_complex operator*(const mp_complex& z, const unsigned short n)     { return mp_complex(z.real() * n, z.imag() * n); }
+      inline mp_complex operator*(const mp_complex& z, const unsigned int n)       { return mp_complex(z.real() * n, z.imag() * n); }
+      inline mp_complex operator*(const mp_complex& z, const unsigned long n)      { return mp_complex(z.real() * n, z.imag() * n); }
+      inline mp_complex operator*(const mp_complex& z, const unsigned long long n) { return mp_complex(z.real() * n, z.imag() * n); }
+      inline mp_complex operator*(const mp_complex& z, const float f)              { return mp_complex(z.real() * f, z.imag() * f); }
+      inline mp_complex operator*(const mp_complex& z, const double d)             { return mp_complex(z.real() * d, z.imag() * d); }
+      inline mp_complex operator*(const mp_complex& z, const long double ld)       { return mp_complex(z.real() * ld, z.imag() * ld); }
+
+      inline mp_complex operator/(const mp_complex& z, const char n)               { return mp_complex(z.real() / n, z.imag() / n); }
+      inline mp_complex operator/(const mp_complex& z, const signed char n)        { return mp_complex(z.real() / n, z.imag() / n); }
+      inline mp_complex operator/(const mp_complex& z, const signed short n)       { return mp_complex(z.real() / n, z.imag() / n); }
+      inline mp_complex operator/(const mp_complex& z, const signed int n)         { return mp_complex(z.real() / n, z.imag() / n); }
+      inline mp_complex operator/(const mp_complex& z, const signed long n)        { return mp_complex(z.real() / n, z.imag() / n); }
+      inline mp_complex operator/(const mp_complex& z, const signed long long n)   { return mp_complex(z.real() / n, z.imag() / n); }
+      inline mp_complex operator/(const mp_complex& z, const unsigned char n)      { return mp_complex(z.real() / n, z.imag() / n); }
+      inline mp_complex operator/(const mp_complex& z, const wchar_t n)            { return mp_complex(z.real() / n, z.imag() / n); }
+      inline mp_complex operator/(const mp_complex& z, const unsigned short n)     { return mp_complex(z.real() / n, z.imag() / n); }
+      inline mp_complex operator/(const mp_complex& z, const unsigned int n)       { return mp_complex(z.real() / n, z.imag() / n); }
+      inline mp_complex operator/(const mp_complex& z, const unsigned long n)      { return mp_complex(z.real() / n, z.imag() / n); }
+      inline mp_complex operator/(const mp_complex& z, const unsigned long long n) { return mp_complex(z.real() / n, z.imag() / n); }
+      inline mp_complex operator/(const mp_complex& z, const float f)              { return mp_complex(z.real() / f, z.imag() / f); }
+      inline mp_complex operator/(const mp_complex& z, const double d)             { return mp_complex(z.real() / d, z.imag() / d); }
+      inline mp_complex operator/(const mp_complex& z, const long double ld)       { return mp_complex(z.real() / ld, z.imag() / ld); }
+
+      // Global add/sub/mul/div of all built-in types with const mp_complex&.
+      inline mp_complex operator+(const char n, const mp_complex& v)               { return mp_complex(v.real() + n, v.imag()); }
+      inline mp_complex operator+(const signed char n, const mp_complex& v)        { return mp_complex(v.real() + n, v.imag()); }
+      inline mp_complex operator+(const signed short n, const mp_complex& v)       { return mp_complex(v.real() + n, v.imag()); }
+      inline mp_complex operator+(const signed int n, const mp_complex& v)         { return mp_complex(v.real() + n, v.imag()); }
+      inline mp_complex operator+(const signed long n, const mp_complex& v)        { return mp_complex(v.real() + n, v.imag()); }
+      inline mp_complex operator+(const signed long long n, const mp_complex& v)   { return mp_complex(v.real() + n, v.imag()); }
+      inline mp_complex operator+(const unsigned char n, const mp_complex& v)      { return mp_complex(v.real() + n, v.imag()); }
+      inline mp_complex operator+(const wchar_t n, const mp_complex& v)            { return mp_complex(v.real() + n, v.imag()); }
+      inline mp_complex operator+(const unsigned short n, const mp_complex& v)     { return mp_complex(v.real() + n, v.imag()); }
+      inline mp_complex operator+(const unsigned int n, const mp_complex& v)       { return mp_complex(v.real() + n, v.imag()); }
+      inline mp_complex operator+(const unsigned long n, const mp_complex& v)      { return mp_complex(v.real() + n, v.imag()); }
+      inline mp_complex operator+(const unsigned long long n, const mp_complex& v) { return mp_complex(v.real() + n, v.imag()); }
+      inline mp_complex operator+(const float f, const mp_complex& v)              { return mp_complex(v.real() + f, v.imag() + f); }
+      inline mp_complex operator+(const double d, const mp_complex& v)             { return mp_complex(v.real() + d, v.imag() + d); }
+      inline mp_complex operator+(const long double ld, const mp_complex& v)       { return mp_complex(v.real() + ld, v.imag() + ld); }
+
+      inline mp_complex operator-(const char n, const mp_complex& v)               { return mp_complex(n - v.real(), -v.imag()); }
+      inline mp_complex operator-(const signed char n, const mp_complex& v)        { return mp_complex(n - v.real(), -v.imag()); }
+      inline mp_complex operator-(const signed short n, const mp_complex& v)       { return mp_complex(n - v.real(), -v.imag()); }
+      inline mp_complex operator-(const signed int n, const mp_complex& v)         { return mp_complex(n - v.real(), -v.imag()); }
+      inline mp_complex operator-(const signed long n, const mp_complex& v)        { return mp_complex(n - v.real(), -v.imag()); }
+      inline mp_complex operator-(const signed long long n, const mp_complex& v)   { return mp_complex(n - v.real(), -v.imag()); }
+      inline mp_complex operator-(const unsigned char n, const mp_complex& v)      { return mp_complex(n - v.real(), -v.imag()); }
+      inline mp_complex operator-(const wchar_t n, const mp_complex& v)            { return mp_complex(n - v.real(), -v.imag()); }
+      inline mp_complex operator-(const unsigned short n, const mp_complex& v)     { return mp_complex(n - v.real(), -v.imag()); }
+      inline mp_complex operator-(const unsigned int n, const mp_complex& v)       { return mp_complex(n - v.real(), -v.imag()); }
+      inline mp_complex operator-(const unsigned long n, const mp_complex& v)      { return mp_complex(n - v.real(), -v.imag()); }
+      inline mp_complex operator-(const unsigned long long n, const mp_complex& v) { return mp_complex(n - v.real(), -v.imag()); }
+      inline mp_complex operator-(const float f, const mp_complex& v)              { return mp_complex(f - v.real(), -v.imag()); }
+      inline mp_complex operator-(const double d, const mp_complex& v)             { return mp_complex(d - v.real(), -v.imag()); }
+      inline mp_complex operator-(const long double ld, const mp_complex& v)       { return mp_complex(ld - v.real(), -v.imag()); }
+
+      inline mp_complex operator*(const char n, const mp_complex& v)               { return mp_complex(v.real() * n, v.imag() * n); }
+      inline mp_complex operator*(const signed char n, const mp_complex& v)        { return mp_complex(v.real() * n, v.imag() * n); }
+      inline mp_complex operator*(const signed short n, const mp_complex& v)       { return mp_complex(v.real() * n, v.imag() * n); }
+      inline mp_complex operator*(const signed int n, const mp_complex& v)         { return mp_complex(v.real() * n, v.imag() * n); }
+      inline mp_complex operator*(const signed long n, const mp_complex& v)        { return mp_complex(v.real() * n, v.imag() * n); }
+      inline mp_complex operator*(const signed long long n, const mp_complex& v)   { return mp_complex(v.real() * n, v.imag() * n); }
+      inline mp_complex operator*(const unsigned char n, const mp_complex& v)      { return mp_complex(v.real() * n, v.imag() * n); }
+      inline mp_complex operator*(const wchar_t n, const mp_complex& v)            { return mp_complex(v.real() * n, v.imag() * n); }
+      inline mp_complex operator*(const unsigned short n, const mp_complex& v)     { return mp_complex(v.real() * n, v.imag() * n); }
+      inline mp_complex operator*(const unsigned int n, const mp_complex& v)       { return mp_complex(v.real() * n, v.imag() * n); }
+      inline mp_complex operator*(const unsigned long n, const mp_complex& v)      { return mp_complex(v.real() * n, v.imag() * n); }
+      inline mp_complex operator*(const unsigned long long n, const mp_complex& v) { return mp_complex(v.real() * n, v.imag() * n); }
+      inline mp_complex operator*(const float f, const mp_complex& v)              { return mp_complex(v.real() * f, v.imag() * f); }
+      inline mp_complex operator*(const double d, const mp_complex& v)             { return mp_complex(v.real() * d, v.imag() * d); }
+      inline mp_complex operator*(const long double ld, const mp_complex& v)       { return mp_complex(v.real() * ld, v.imag() * ld); }
+
+      inline mp_complex operator/(const char n, const mp_complex& v)               { const mp_float v_norm = v.norm(); return mp_complex((n * v.real()) / v_norm, -(n * v.imag()) / v_norm); }
+      inline mp_complex operator/(const signed char n, const mp_complex& v)        { const mp_float v_norm = v.norm(); return mp_complex((n * v.real()) / v_norm, -(n * v.imag()) / v_norm); }
+      inline mp_complex operator/(const signed short n, const mp_complex& v)       { const mp_float v_norm = v.norm(); return mp_complex((n * v.real()) / v_norm, -(n * v.imag()) / v_norm); }
+      inline mp_complex operator/(const signed int n, const mp_complex& v)         { const mp_float v_norm = v.norm(); return mp_complex((n * v.real()) / v_norm, -(n * v.imag()) / v_norm); }
+      inline mp_complex operator/(const signed long n, const mp_complex& v)        { const mp_float v_norm = v.norm(); return mp_complex((n * v.real()) / v_norm, -(n * v.imag()) / v_norm); }
+      inline mp_complex operator/(const signed long long n, const mp_complex& v)   { const mp_float v_norm = v.norm(); return mp_complex((n * v.real()) / v_norm, -(n * v.imag()) / v_norm); }
+      inline mp_complex operator/(const unsigned char n, const mp_complex& v)      { const mp_float v_norm = v.norm(); return mp_complex((n * v.real()) / v_norm, -(n * v.imag()) / v_norm); }
+      inline mp_complex operator/(const wchar_t n, const mp_complex& v)            { const mp_float v_norm = v.norm(); return mp_complex((n * v.real()) / v_norm, -(n * v.imag()) / v_norm); }
+      inline mp_complex operator/(const unsigned short n, const mp_complex& v)     { const mp_float v_norm = v.norm(); return mp_complex((n * v.real()) / v_norm, -(n * v.imag()) / v_norm); }
+      inline mp_complex operator/(const unsigned int n, const mp_complex& v)       { const mp_float v_norm = v.norm(); return mp_complex((n * v.real()) / v_norm, -(n * v.imag()) / v_norm); }
+      inline mp_complex operator/(const unsigned long n, const mp_complex& v)      { const mp_float v_norm = v.norm(); return mp_complex((n * v.real()) / v_norm, -(n * v.imag()) / v_norm); }
+      inline mp_complex operator/(const unsigned long long n, const mp_complex& v) { const mp_float v_norm = v.norm(); return mp_complex((n * v.real()) / v_norm, -(n * v.imag()) / v_norm); }
+      inline mp_complex operator/(const float f, const mp_complex& v)              { const mp_float v_norm = v.norm(); return mp_complex((f * v.real()) / v_norm, -(f * v.imag()) / v_norm); }
+      inline mp_complex operator/(const double d, const mp_complex& v)             { const mp_float v_norm = v.norm(); return mp_complex((d * v.real()) / v_norm, -(d * v.imag()) / v_norm); }
+      inline mp_complex operator/(const long double ld, const mp_complex& v)       { const mp_float v_norm = v.norm(); return mp_complex((ld * v.real()) / v_norm, -(ld * v.imag()) / v_norm); }
+
+      // Global self add/sub/mul/div of mp_complex& with all built-in types.
+      inline mp_complex& operator+=(mp_complex& z, const char n)               { z.Re += n; return z; }
+      inline mp_complex& operator+=(mp_complex& z, const signed char n)        { z.Re += n; return z; }
+      inline mp_complex& operator+=(mp_complex& z, const signed short n)       { z.Re += n; return z; }
+      inline mp_complex& operator+=(mp_complex& z, const signed int n)         { z.Re += n; return z; }
+      inline mp_complex& operator+=(mp_complex& z, const signed long n)        { z.Re += n; return z; }
+      inline mp_complex& operator+=(mp_complex& z, const signed long long n)   { z.Re += n; return z; }
+      inline mp_complex& operator+=(mp_complex& z, const unsigned char n)      { z.Re += n; return z; }
+      inline mp_complex& operator+=(mp_complex& z, const wchar_t n)            { z.Re += n; return z; }
+      inline mp_complex& operator+=(mp_complex& z, const unsigned short n)     { z.Re += n; return z; }
+      inline mp_complex& operator+=(mp_complex& z, const unsigned int n)       { z.Re += n; return z; }
+      inline mp_complex& operator+=(mp_complex& z, const unsigned long n)      { z.Re += n; return z; }
+      inline mp_complex& operator+=(mp_complex& z, const unsigned long long n) { z.Re += n; return z; }
+      inline mp_complex& operator+=(mp_complex& z, const float f)              { z.Re += f; return z; }
+      inline mp_complex& operator+=(mp_complex& z, const double d)             { z.Re += d; return z; }
+      inline mp_complex& operator+=(mp_complex& z, const long double ld)       { z.Re += ld; return z; }
+
+      inline mp_complex& operator-=(mp_complex& z, const char n)               { z.Re -= n; return z; }
+      inline mp_complex& operator-=(mp_complex& z, const signed char n)        { z.Re -= n; return z; }
+      inline mp_complex& operator-=(mp_complex& z, const signed short n)       { z.Re -= n; return z; }
+      inline mp_complex& operator-=(mp_complex& z, const signed int n)         { z.Re -= n; return z; }
+      inline mp_complex& operator-=(mp_complex& z, const signed long n)        { z.Re -= n; return z; }
+      inline mp_complex& operator-=(mp_complex& z, const signed long long n)   { z.Re -= n; return z; }
+      inline mp_complex& operator-=(mp_complex& z, const unsigned char n)      { z.Re -= n; return z; }
+      inline mp_complex& operator-=(mp_complex& z, const wchar_t n)            { z.Re -= n; return z; }
+      inline mp_complex& operator-=(mp_complex& z, const unsigned short n)     { z.Re -= n; return z; }
+      inline mp_complex& operator-=(mp_complex& z, const unsigned int n)       { z.Re -= n; return z; }
+      inline mp_complex& operator-=(mp_complex& z, const unsigned long n)      { z.Re -= n; return z; }
+      inline mp_complex& operator-=(mp_complex& z, const unsigned long long n) { z.Re -= n; return z; }
+      inline mp_complex& operator-=(mp_complex& z, const float f)              { z.Re -= f; return z; }
+      inline mp_complex& operator-=(mp_complex& z, const double d)             { z.Re -= d; return z; }
+      inline mp_complex& operator-=(mp_complex& z, const long double ld)       { z.Re -= ld; return z; }
+
+      inline mp_complex& operator*=(mp_complex& z, const char n)               { z.Re *= n; z.Im *= n; return z; }
+      inline mp_complex& operator*=(mp_complex& z, const signed char n)        { z.Re *= n; z.Im *= n; return z; }
+      inline mp_complex& operator*=(mp_complex& z, const signed short n)       { z.Re *= n; z.Im *= n; return z; }
+      inline mp_complex& operator*=(mp_complex& z, const signed int n)         { z.Re *= n; z.Im *= n; return z; }
+      inline mp_complex& operator*=(mp_complex& z, const signed long n)        { z.Re *= n; z.Im *= n; return z; }
+      inline mp_complex& operator*=(mp_complex& z, const signed long long n)   { z.Re *= n; z.Im *= n; return z; }
+      inline mp_complex& operator*=(mp_complex& z, const unsigned char n)      { z.Re *= n; z.Im *= n; return z; }
+      inline mp_complex& operator*=(mp_complex& z, const wchar_t n)            { z.Re *= n; z.Im *= n; return z; }
+      inline mp_complex& operator*=(mp_complex& z, const unsigned short n)     { z.Re *= n; z.Im *= n; return z; }
+      inline mp_complex& operator*=(mp_complex& z, const unsigned int n)       { z.Re *= n; z.Im *= n; return z; }
+      inline mp_complex& operator*=(mp_complex& z, const unsigned long n)      { z.Re *= n; z.Im *= n; return z; }
+      inline mp_complex& operator*=(mp_complex& z, const unsigned long long n) { z.Re *= n; z.Im *= n; return z; }
+      inline mp_complex& operator*=(mp_complex& z, const float f)              { z.Re *= f; z.Im *= f; return z; }
+      inline mp_complex& operator*=(mp_complex& z, const double d)             { z.Re *= d; z.Im *= d; return z; }
+      inline mp_complex& operator*=(mp_complex& z, const long double ld)       { z.Re *= ld; z.Im *= ld; return z; }
+
+      inline mp_complex& operator/=(mp_complex& z, const char n)               { z.Re /= n; z.Im /= n; return z; }
+      inline mp_complex& operator/=(mp_complex& z, const signed char n)        { z.Re /= n; z.Im /= n; return z; }
+      inline mp_complex& operator/=(mp_complex& z, const signed short n)       { z.Re /= n; z.Im /= n; return z; }
+      inline mp_complex& operator/=(mp_complex& z, const signed int n)         { z.Re /= n; z.Im /= n; return z; }
+      inline mp_complex& operator/=(mp_complex& z, const signed long n)        { z.Re /= n; z.Im /= n; return z; }
+      inline mp_complex& operator/=(mp_complex& z, const signed long long n)   { z.Re /= n; z.Im /= n; return z; }
+      inline mp_complex& operator/=(mp_complex& z, const unsigned char n)      { z.Re /= n; z.Im /= n; return z; }
+      inline mp_complex& operator/=(mp_complex& z, const wchar_t n)            { z.Re /= n; z.Im /= n; return z; }
+      inline mp_complex& operator/=(mp_complex& z, const unsigned short n)     { z.Re /= n; z.Im /= n; return z; }
+      inline mp_complex& operator/=(mp_complex& z, const unsigned int n)       { z.Re /= n; z.Im /= n; return z; }
+      inline mp_complex& operator/=(mp_complex& z, const unsigned long n)      { z.Re /= n; z.Im /= n; return z; }
+      inline mp_complex& operator/=(mp_complex& z, const unsigned long long n) { z.Re /= n; z.Im /= n; return z; }
+      inline mp_complex& operator/=(mp_complex& z, const float f)              { z.Re /= f; z.Im /= f; return z; }
+      inline mp_complex& operator/=(mp_complex& z, const double d)             { z.Re /= d; z.Im /= d; return z; }
+      inline mp_complex& operator/=(mp_complex& z, const long double ld)       { z.Re /= ld; z.Im /= ld; return z; }
+
+      // Global self add/sub/mul/div of mp_complex& with const mp_float&.
+      inline mp_complex& operator+=(mp_complex& z, const mp_float& v) { z.Re += v;            return z; }
+      inline mp_complex& operator-=(mp_complex& z, const mp_float& v) { z.Re -= v;            return z; }
+      inline mp_complex& operator*=(mp_complex& z, const mp_float& v) { z.Re *= v; z.Im *= v; return z; }
+      inline mp_complex& operator/=(mp_complex& z, const mp_float& v) { z.Re /= v; z.Im /= v; return z; }
+    }
+  }
+
+#endif // _MP_COMPLEX_HPP_

Added: sandbox/multiprecision/boost/mp_float.hpp
==============================================================================
--- (empty file)
+++ sandbox/multiprecision/boost/mp_float.hpp	2011-09-21 17:00:15 EDT (Wed, 21 Sep 2011)
@@ -0,0 +1,522 @@
+
+//          Copyright Christopher Kormanyos 2002 - 2011.
+// Distributed under the Boost Software License, Version 1.0.
+//    (See accompanying file LICENSE_1_0.txt or copy at
+//          http://www.boost.org/LICENSE_1_0.txt)
+
+// This work is based on an earlier work:
+// "Algorithm 910: A Portable C++ Multiple-Precision System for Special-Function Calculations",
+// in ACM TOMS, {VOL 37, ISSUE 4, (February 2011)} (C) ACM, 2011. http://doi.acm.org/10.1145/1916461.1916469
+
+#ifndef _MP_FLOAT_2004_02_08_HPP_
+  #define _MP_FLOAT_2004_02_08_HPP_
+
+  #include <limits>
+
+  // Select the mp_float back-end big-number type
+  // by defining BOOST_MULTIPRECISION_BACKEND_MP_FLOAT_TYPE_xxx, for example,
+  // as a compiler defined preprocessor macro.
+
+  #if defined(BOOST_MULTIPRECISION_BACKEND_MP_FLOAT_TYPE_EFX)
+    #include <boost/multiprecision/mp_float_efx.hpp>
+  #elif defined(BOOST_MULTIPRECISION_BACKEND_MP_FLOAT_TYPE_GMP)
+    #include <boost/multiprecision/mp_float_gmp.hpp>
+  #elif defined(BOOST_MULTIPRECISION_BACKEND_MP_FLOAT_TYPE_MPFR)
+    #include <boost/multiprecision/mp_float_mpfr.hpp>
+  #else
+    #error mp_float type undefined!
+  #endif
+
+  namespace boost
+  {
+    namespace multiprecision
+    {
+      const mp_float& zero     (void);
+      const mp_float& one      (void);
+      const mp_float& half     (void);
+      const mp_float& value_min(void);
+      const mp_float& value_max(void);
+      const mp_float& value_eps(void);
+      const mp_float& value_inf(void);
+      const mp_float& value_nan(void);
+    }
+  }
+
+  // Specialization of std::numeric_limits<mp_float>.
+  namespace std
+  {
+    template <> class numeric_limits<boost::multiprecision::mp_float>
+    {
+    public:
+      static const bool                    is_specialized    = true;
+      static const bool                    is_signed         = true;
+      static const bool                    is_integer        = false;
+      static const bool                    is_exact          = false;
+      static const bool                    is_bounded        = true;
+      static const bool                    is_modulo         = false;
+      static const bool                    is_iec559         = false;
+      static const int                     digits            = boost::multiprecision::mp_float::mp_digits;
+      static const int                     digits10          = boost::multiprecision::mp_float::mp_digits10;
+      static const int                     max_digits10      = boost::multiprecision::mp_float::mp_max_digits10;
+      static const boost::int64_t                   min_exponent      = boost::multiprecision::mp_float::mp_min_exp;      // Type differs from int.
+      static const boost::int64_t                   min_exponent10    = boost::multiprecision::mp_float::mp_min_exp10;    // Type differs from int.
+      static const boost::int64_t                   max_exponent      = boost::multiprecision::mp_float::mp_max_exp;      // Type differs from int.
+      static const boost::int64_t                   max_exponent10    = boost::multiprecision::mp_float::mp_max_exp10;    // Type differs from int.
+      static const int                     radix             = boost::multiprecision::mp_float::mp_radix;
+      static const std::float_round_style  round_style       = std::round_to_nearest;
+      static const bool                    has_infinity      = true;
+      static const bool                    has_quiet_NaN     = true;
+      static const bool                    has_signaling_NaN = false;
+      static const std::float_denorm_style has_denorm        = std::denorm_absent;
+      static const bool                    has_denorm_loss   = false;
+      static const bool                    traps             = false;
+      static const bool                    tinyness_before   = false;
+
+      static const boost::multiprecision::mp_float& (min)        (void) throw() { return boost::multiprecision::value_min(); }
+      static const boost::multiprecision::mp_float& (max)        (void) throw() { return boost::multiprecision::value_max(); }
+      static const boost::multiprecision::mp_float& lowest       (void) throw() { return boost::multiprecision::zero(); }
+      static const boost::multiprecision::mp_float& epsilon      (void) throw() { return boost::multiprecision::value_eps(); }
+      static const boost::multiprecision::mp_float& round_error  (void) throw() { return boost::multiprecision::half(); }
+      static const boost::multiprecision::mp_float& infinity     (void) throw() { return boost::multiprecision::value_inf(); }
+      static const boost::multiprecision::mp_float& quiet_NaN    (void) throw() { return boost::multiprecision::value_nan(); }
+      static const boost::multiprecision::mp_float& signaling_NaN(void) throw() { return boost::multiprecision::zero(); }
+      static const boost::multiprecision::mp_float& denorm_min   (void) throw() { return boost::multiprecision::zero(); }
+    };
+  }
+
+  namespace boost
+  {
+    namespace multiprecision
+    {
+      // Global operators post-increment and post-decrement.
+      inline mp_float operator++(mp_float& u, int) { const mp_float v(u); ++u; return v; }
+      inline mp_float operator--(mp_float& u, int) { const mp_float v(u); --u; return v; }
+
+      // Global unary operators of mp_float reference.
+      inline       mp_float  operator-(const mp_float& u) { return mp_float(u).negate(); }
+      inline       mp_float& operator+(      mp_float& u) { return u; }
+      inline const mp_float& operator+(const mp_float& u) { return u; }
+
+      // Global add/sub/mul/div of const mp_float& with const mp_float&.
+      inline mp_float operator+(const mp_float& u, const mp_float& v) { return mp_float(u) += v; }
+      inline mp_float operator-(const mp_float& u, const mp_float& v) { return mp_float(u) -= v; }
+      inline mp_float operator*(const mp_float& u, const mp_float& v) { return mp_float(u) *= v; }
+      inline mp_float operator/(const mp_float& u, const mp_float& v) { return mp_float(u) /= v; }
+
+      // Global add/sub/mul/div of const mp_float& with all built-in types.
+      inline mp_float operator+(const mp_float& u, const char n)               { return ((!std::numeric_limits<char>::is_signed) ? mp_float(u).add_unsigned_long_long(static_cast<unsigned long long>(n)) : mp_float(u).add_signed_long_long(static_cast<signed long long>(n))); }
+      inline mp_float operator+(const mp_float& u, const signed char n)        { return mp_float(u).add_signed_long_long(static_cast<signed long long>(n)); }
+      inline mp_float operator+(const mp_float& u, const signed short n)       { return mp_float(u).add_signed_long_long(static_cast<signed long long>(n)); }
+      inline mp_float operator+(const mp_float& u, const signed int n)         { return mp_float(u).add_signed_long_long(static_cast<signed long long>(n)); }
+      inline mp_float operator+(const mp_float& u, const signed long n)        { return mp_float(u).add_signed_long_long(static_cast<signed long long>(n)); }
+      inline mp_float operator+(const mp_float& u, const signed long long n)   { return mp_float(u).add_signed_long_long(static_cast<signed long long>(n)); }
+      inline mp_float operator+(const mp_float& u, const unsigned char n)      { return mp_float(u).add_unsigned_long_long(static_cast<unsigned long long>(n)); }
+      inline mp_float operator+(const mp_float& u, const wchar_t n)            { return ((!std::numeric_limits<wchar_t>::is_signed) ? mp_float(u).add_unsigned_long_long(static_cast<unsigned long long>(n)) : mp_float(u).add_signed_long_long(static_cast<signed long long>(n))); }
+      inline mp_float operator+(const mp_float& u, const unsigned short n)     { return mp_float(u).add_unsigned_long_long(static_cast<unsigned long long>(n)); }
+      inline mp_float operator+(const mp_float& u, const unsigned int n)       { return mp_float(u).add_unsigned_long_long(static_cast<unsigned long long>(n)); }
+      inline mp_float operator+(const mp_float& u, const unsigned long n)      { return mp_float(u).add_unsigned_long_long(static_cast<unsigned long long>(n)); }
+      inline mp_float operator+(const mp_float& u, const unsigned long long n) { return mp_float(u).add_unsigned_long_long(static_cast<unsigned long long>(n)); }
+      inline mp_float operator+(const mp_float& u, const float f)              { return mp_float(u) += mp_float(f); }
+      inline mp_float operator+(const mp_float& u, const double d)             { return mp_float(u) += mp_float(d); }
+      inline mp_float operator+(const mp_float& u, const long double ld)       { return mp_float(u) += mp_float(ld); }
+
+      inline mp_float operator-(const mp_float& u, const char n)               { return ((!std::numeric_limits<char>::is_signed) ? mp_float(u).sub_unsigned_long_long(static_cast<unsigned long long>(n)) : mp_float(u).sub_signed_long_long(static_cast<signed long long>(n))); }
+      inline mp_float operator-(const mp_float& u, const signed char n)        { return mp_float(u).sub_signed_long_long(static_cast<signed long long>(n)); }
+      inline mp_float operator-(const mp_float& u, const signed short n)       { return mp_float(u).sub_signed_long_long(static_cast<signed long long>(n)); }
+      inline mp_float operator-(const mp_float& u, const signed int n)         { return mp_float(u).sub_signed_long_long(static_cast<signed long long>(n)); }
+      inline mp_float operator-(const mp_float& u, const signed long n)        { return mp_float(u).sub_signed_long_long(static_cast<signed long long>(n)); }
+      inline mp_float operator-(const mp_float& u, const signed long long n)   { return mp_float(u).sub_signed_long_long(static_cast<signed long long>(n)); }
+      inline mp_float operator-(const mp_float& u, const unsigned char n)      { return mp_float(u).sub_unsigned_long_long(static_cast<unsigned long long>(n)); }
+      inline mp_float operator-(const mp_float& u, const wchar_t n)            { return ((!std::numeric_limits<wchar_t>::is_signed) ? mp_float(u).sub_unsigned_long_long(static_cast<unsigned long long>(n)) : mp_float(u).sub_signed_long_long(static_cast<signed long long>(n))); }
+      inline mp_float operator-(const mp_float& u, const unsigned short n)     { return mp_float(u).sub_unsigned_long_long(static_cast<unsigned long long>(n)); }
+      inline mp_float operator-(const mp_float& u, const unsigned int n)       { return mp_float(u).sub_unsigned_long_long(static_cast<unsigned long long>(n)); }
+      inline mp_float operator-(const mp_float& u, const unsigned long n)      { return mp_float(u).sub_unsigned_long_long(static_cast<unsigned long long>(n)); }
+      inline mp_float operator-(const mp_float& u, const unsigned long long n) { return mp_float(u).sub_unsigned_long_long(static_cast<unsigned long long>(n)); }
+      inline mp_float operator-(const mp_float& u, const float f)              { return mp_float(u) -= mp_float(f); }
+      inline mp_float operator-(const mp_float& u, const double d)             { return mp_float(u) -= mp_float(d); }
+      inline mp_float operator-(const mp_float& u, const long double ld)       { return mp_float(u) -= mp_float(ld); }
+
+      inline mp_float operator*(const mp_float& u, const char n)               { return ((!std::numeric_limits<char>::is_signed) ? mp_float(u).mul_unsigned_long_long(static_cast<unsigned long long>(n)) : mp_float(u).mul_signed_long_long(static_cast<signed long long>(n))); }
+      inline mp_float operator*(const mp_float& u, const signed char n)        { return mp_float(u).mul_signed_long_long(static_cast<signed long long>(n)); }
+      inline mp_float operator*(const mp_float& u, const signed short n)       { return mp_float(u).mul_signed_long_long(static_cast<signed long long>(n)); }
+      inline mp_float operator*(const mp_float& u, const signed int n)         { return mp_float(u).mul_signed_long_long(static_cast<signed long long>(n)); }
+      inline mp_float operator*(const mp_float& u, const signed long n)        { return mp_float(u).mul_signed_long_long(static_cast<signed long long>(n)); }
+      inline mp_float operator*(const mp_float& u, const signed long long n)   { return mp_float(u).mul_signed_long_long(static_cast<signed long long>(n)); }
+      inline mp_float operator*(const mp_float& u, const unsigned char n)      { return mp_float(u).mul_unsigned_long_long(static_cast<unsigned long long>(n)); }
+      inline mp_float operator*(const mp_float& u, const wchar_t n)            { return ((!std::numeric_limits<wchar_t>::is_signed) ? mp_float(u).mul_unsigned_long_long(static_cast<unsigned long long>(n)) : mp_float(u).mul_signed_long_long(static_cast<signed long long>(n))); }
+      inline mp_float operator*(const mp_float& u, const unsigned short n)     { return mp_float(u).mul_unsigned_long_long(static_cast<unsigned long long>(n)); }
+      inline mp_float operator*(const mp_float& u, const unsigned int n)       { return mp_float(u).mul_unsigned_long_long(static_cast<unsigned long long>(n)); }
+      inline mp_float operator*(const mp_float& u, const unsigned long n)      { return mp_float(u).mul_unsigned_long_long(static_cast<unsigned long long>(n)); }
+      inline mp_float operator*(const mp_float& u, const unsigned long long n) { return mp_float(u).mul_unsigned_long_long(static_cast<unsigned long long>(n)); }
+      inline mp_float operator*(const mp_float& u, const float f)              { return mp_float(u) *= mp_float(f); }
+      inline mp_float operator*(const mp_float& u, const double d)             { return mp_float(u) *= mp_float(d); }
+      inline mp_float operator*(const mp_float& u, const long double ld)       { return mp_float(u) *= mp_float(ld); }
+
+      inline mp_float operator/(const mp_float& u, const char n)               { return ((!std::numeric_limits<char>::is_signed) ? mp_float(u).div_unsigned_long_long(static_cast<unsigned long long>(n)) : mp_float(u).div_signed_long_long(static_cast<signed long long>(n))); }
+      inline mp_float operator/(const mp_float& u, const signed char n)        { return mp_float(u).div_signed_long_long(static_cast<signed long long>(n)); }
+      inline mp_float operator/(const mp_float& u, const signed short n)       { return mp_float(u).div_signed_long_long(static_cast<signed long long>(n)); }
+      inline mp_float operator/(const mp_float& u, const signed int n)         { return mp_float(u).div_signed_long_long(static_cast<signed long long>(n)); }
+      inline mp_float operator/(const mp_float& u, const signed long n)        { return mp_float(u).div_signed_long_long(static_cast<signed long long>(n)); }
+      inline mp_float operator/(const mp_float& u, const signed long long n)   { return mp_float(u).div_signed_long_long(static_cast<signed long long>(n)); }
+      inline mp_float operator/(const mp_float& u, const unsigned char n)      { return mp_float(u).div_unsigned_long_long(static_cast<unsigned long long>(n)); }
+      inline mp_float operator/(const mp_float& u, const wchar_t n)            { return ((!std::numeric_limits<wchar_t>::is_signed) ? mp_float(u).div_unsigned_long_long(static_cast<unsigned long long>(n)) : mp_float(u).div_signed_long_long(static_cast<signed long long>(n))); }
+      inline mp_float operator/(const mp_float& u, const unsigned short n)     { return mp_float(u).div_unsigned_long_long(static_cast<unsigned long long>(n)); }
+      inline mp_float operator/(const mp_float& u, const unsigned int n)       { return mp_float(u).div_unsigned_long_long(static_cast<unsigned long long>(n)); }
+      inline mp_float operator/(const mp_float& u, const unsigned long n)      { return mp_float(u).div_unsigned_long_long(static_cast<unsigned long long>(n)); }
+      inline mp_float operator/(const mp_float& u, const unsigned long long n) { return mp_float(u).div_unsigned_long_long(static_cast<unsigned long long>(n)); }
+      inline mp_float operator/(const mp_float& u, const float f)              { return mp_float(u) /= mp_float(f); }
+      inline mp_float operator/(const mp_float& u, const double d)             { return mp_float(u) /= mp_float(d); }
+      inline mp_float operator/(const mp_float& u, const long double ld)       { return mp_float(u) /= mp_float(ld); }
+
+      // Global add/sub/mul/div of all built-in types with const mp_float&.
+      inline mp_float operator+(const char n, const mp_float& u)               { return ((!std::numeric_limits<char>::is_signed) ? mp_float(u).add_unsigned_long_long(n) : mp_float(u).add_signed_long_long(n)); }
+      inline mp_float operator+(const signed char n, const mp_float& u)        { return mp_float(u).add_signed_long_long(n); }
+      inline mp_float operator+(const signed short n, const mp_float& u)       { return mp_float(u).add_signed_long_long(n); }
+      inline mp_float operator+(const signed int n, const mp_float& u)         { return mp_float(u).add_signed_long_long(n); }
+      inline mp_float operator+(const signed long n, const mp_float& u)        { return mp_float(u).add_signed_long_long(n); }
+      inline mp_float operator+(const signed long long n, const mp_float& u)   { return mp_float(u).add_signed_long_long(n); }
+      inline mp_float operator+(const unsigned char n, const mp_float& u)      { return mp_float(u).add_unsigned_long_long(n); }
+      inline mp_float operator+(const wchar_t n, const mp_float& u)            { return ((!std::numeric_limits<wchar_t>::is_signed) ? mp_float(u).add_unsigned_long_long(n) : mp_float(u).add_signed_long_long(n)); }
+      inline mp_float operator+(const unsigned short n, const mp_float& u)     { return mp_float(u).add_unsigned_long_long(n); }
+      inline mp_float operator+(const unsigned int n, const mp_float& u)       { return mp_float(u).add_unsigned_long_long(n); }
+      inline mp_float operator+(const unsigned long n, const mp_float& u)      { return mp_float(u).add_unsigned_long_long(n); }
+      inline mp_float operator+(const unsigned long long n, const mp_float& u) { return mp_float(u).add_unsigned_long_long(n); }
+      inline mp_float operator+(const float f, const mp_float& u)              { return mp_float(f) += u; }
+      inline mp_float operator+(const double d, const mp_float& u)             { return mp_float(d) += u; }
+      inline mp_float operator+(const long double ld, const mp_float& u)       { return mp_float(ld) += u; }
+
+      inline mp_float operator-(const char n, const mp_float& u)               { return mp_float(n)  -= u; }
+      inline mp_float operator-(const signed char n, const mp_float& u)        { return mp_float(n)  -= u; }
+      inline mp_float operator-(const signed short n, const mp_float& u)       { return mp_float(n)  -= u; }
+      inline mp_float operator-(const signed int n, const mp_float& u)         { return mp_float(n)  -= u; }
+      inline mp_float operator-(const signed long n, const mp_float& u)        { return mp_float(n)  -= u; }
+      inline mp_float operator-(const signed long long n, const mp_float& u)   { return mp_float(n)  -= u; }
+      inline mp_float operator-(const unsigned char n, const mp_float& u)      { return mp_float(n)  -= u; }
+      inline mp_float operator-(const wchar_t n, const mp_float& u)            { return mp_float(n)  -= u; }
+      inline mp_float operator-(const unsigned short n, const mp_float& u)     { return mp_float(n)  -= u; }
+      inline mp_float operator-(const unsigned int n, const mp_float& u)       { return mp_float(n)  -= u; }
+      inline mp_float operator-(const unsigned long n, const mp_float& u)      { return mp_float(n)  -= u; }
+      inline mp_float operator-(const unsigned long long n, const mp_float& u) { return mp_float(n)  -= u; }
+      inline mp_float operator-(const float f, const mp_float& u)              { return mp_float(f)  -= u; }
+      inline mp_float operator-(const double d, const mp_float& u)             { return mp_float(d)  -= u; }
+      inline mp_float operator-(const long double ld, const mp_float& u)       { return mp_float(ld) -= u; }
+
+      inline mp_float operator*(const char n, const mp_float& u)               { return ((!std::numeric_limits<char>::is_signed) ? mp_float(u).mul_unsigned_long_long(n) : mp_float(u).mul_signed_long_long(n)); }
+      inline mp_float operator*(const signed char n, const mp_float& u)        { return mp_float(u).mul_signed_long_long(n); }
+      inline mp_float operator*(const signed short n, const mp_float& u)       { return mp_float(u).mul_signed_long_long(n); }
+      inline mp_float operator*(const signed int n, const mp_float& u)         { return mp_float(u).mul_signed_long_long(n); }
+      inline mp_float operator*(const signed long n, const mp_float& u)        { return mp_float(u).mul_signed_long_long(n); }
+      inline mp_float operator*(const signed long long n, const mp_float& u)   { return mp_float(u).mul_signed_long_long(n); }
+      inline mp_float operator*(const unsigned char n, const mp_float& u)      { return mp_float(u).mul_unsigned_long_long(n); }
+      inline mp_float operator*(const wchar_t n, const mp_float& u)            { return ((!std::numeric_limits<wchar_t>::is_signed) ? mp_float(u).mul_unsigned_long_long(n) : mp_float(u).mul_signed_long_long(n)); }
+      inline mp_float operator*(const unsigned short n, const mp_float& u)     { return mp_float(u).mul_unsigned_long_long(n); }
+      inline mp_float operator*(const unsigned int n, const mp_float& u)       { return mp_float(u).mul_unsigned_long_long(n); }
+      inline mp_float operator*(const unsigned long n, const mp_float& u)      { return mp_float(u).mul_unsigned_long_long(n); }
+      inline mp_float operator*(const unsigned long long n, const mp_float& u) { return mp_float(u).mul_unsigned_long_long(n); }
+      inline mp_float operator*(const float f, const mp_float& u)              { return mp_float(f) *= u; }
+      inline mp_float operator*(const double d, const mp_float& u)             { return mp_float(d) *= u; }
+      inline mp_float operator*(const long double ld, const mp_float& u)       { return mp_float(ld) *= u; }
+
+      inline mp_float operator/(const char n, const mp_float& u)               { return mp_float(n)  /= u; }
+      inline mp_float operator/(const signed char n, const mp_float& u)        { return mp_float(n)  /= u; }
+      inline mp_float operator/(const signed short n, const mp_float& u)       { return mp_float(n)  /= u; }
+      inline mp_float operator/(const signed int n, const mp_float& u)         { return mp_float(n)  /= u; }
+      inline mp_float operator/(const signed long n, const mp_float& u)        { return mp_float(n)  /= u; }
+      inline mp_float operator/(const signed long long n, const mp_float& u)   { return mp_float(n)  /= u; }
+      inline mp_float operator/(const unsigned char n, const mp_float& u)      { return mp_float(n)  /= u; }
+      inline mp_float operator/(const wchar_t n, const mp_float& u)            { return mp_float(n)  /= u; }
+      inline mp_float operator/(const unsigned short n, const mp_float& u)     { return mp_float(n)  /= u; }
+      inline mp_float operator/(const unsigned int n, const mp_float& u)       { return mp_float(n)  /= u; }
+      inline mp_float operator/(const unsigned long n, const mp_float& u)      { return mp_float(n)  /= u; }
+      inline mp_float operator/(const unsigned long long n, const mp_float& u) { return mp_float(n)  /= u; }
+      inline mp_float operator/(const float f, const mp_float& u)              { return mp_float(f)  /= u; }
+      inline mp_float operator/(const double d, const mp_float& u)             { return mp_float(d)  /= u; }
+      inline mp_float operator/(const long double ld, const mp_float& u)       { return mp_float(ld) /= u; }
+
+      // Global self add/sub/mul/div of mp_float& with all built-in types.
+      inline mp_float& operator+=(mp_float& u, const char n)               { return ((!std::numeric_limits<char>::is_signed) ? u.add_unsigned_long_long(static_cast<unsigned long long>(n)) : u.add_signed_long_long(static_cast<signed long long>(n))); }
+      inline mp_float& operator+=(mp_float& u, const signed char n)        { return u.add_signed_long_long(n); }
+      inline mp_float& operator+=(mp_float& u, const signed short n)       { return u.add_signed_long_long(n); }
+      inline mp_float& operator+=(mp_float& u, const signed int n)         { return u.add_signed_long_long(n); }
+      inline mp_float& operator+=(mp_float& u, const signed long n)        { return u.add_signed_long_long(n); }
+      inline mp_float& operator+=(mp_float& u, const signed long long n)   { return u.add_signed_long_long(n); }
+      inline mp_float& operator+=(mp_float& u, const unsigned char n)      { return u.add_unsigned_long_long(n); }
+      inline mp_float& operator+=(mp_float& u, const wchar_t n)            { return ((!std::numeric_limits<wchar_t>::is_signed) ? u.add_unsigned_long_long(static_cast<unsigned long long>(n)) : u.add_signed_long_long(static_cast<signed long long>(n))); }
+      inline mp_float& operator+=(mp_float& u, const unsigned short n)     { return u.add_unsigned_long_long(n); }
+      inline mp_float& operator+=(mp_float& u, const unsigned int n)       { return u.add_unsigned_long_long(n); }
+      inline mp_float& operator+=(mp_float& u, const unsigned long n)      { return u.add_unsigned_long_long(n); }
+      inline mp_float& operator+=(mp_float& u, const unsigned long long n) { return u.add_unsigned_long_long(n); }
+      inline mp_float& operator+=(mp_float& u, const float f)              { return u += mp_float(f); }
+      inline mp_float& operator+=(mp_float& u, const double d)             { return u += mp_float(d); }
+      inline mp_float& operator+=(mp_float& u, const long double ld)       { return u += mp_float(ld); }
+
+      inline mp_float& operator-=(mp_float& u, const signed char n)        { return ((!std::numeric_limits<char>::is_signed) ? u.sub_unsigned_long_long(static_cast<unsigned long long>(n)) : u.sub_signed_long_long(static_cast<signed long long>(n))); }
+      inline mp_float& operator-=(mp_float& u, const signed short n)       { return u.sub_signed_long_long(n); }
+      inline mp_float& operator-=(mp_float& u, const signed int n)         { return u.sub_signed_long_long(n); }
+      inline mp_float& operator-=(mp_float& u, const signed long n)        { return u.sub_signed_long_long(n); }
+      inline mp_float& operator-=(mp_float& u, const signed long long n)   { return u.sub_signed_long_long(n); }
+      inline mp_float& operator-=(mp_float& u, const unsigned char n)      { return u.sub_unsigned_long_long(n); }
+      inline mp_float& operator-=(mp_float& u, const wchar_t n)            { return ((!std::numeric_limits<wchar_t>::is_signed) ? u.sub_unsigned_long_long(static_cast<unsigned long long>(n)) : u.sub_signed_long_long(static_cast<signed long long>(n))); }
+      inline mp_float& operator-=(mp_float& u, const unsigned short n)     { return u.sub_unsigned_long_long(n); }
+      inline mp_float& operator-=(mp_float& u, const unsigned int n)       { return u.sub_unsigned_long_long(n); }
+      inline mp_float& operator-=(mp_float& u, const unsigned long n)      { return u.sub_unsigned_long_long(n); }
+      inline mp_float& operator-=(mp_float& u, const unsigned long long n) { return u.sub_unsigned_long_long(n); }
+      inline mp_float& operator-=(mp_float& u, const float f)              { return u -= mp_float(f); }
+      inline mp_float& operator-=(mp_float& u, const double d)             { return u -= mp_float(d); }
+      inline mp_float& operator-=(mp_float& u, const long double ld)       { return u -= mp_float(ld); }
+
+      inline mp_float& operator*=(mp_float& u, const char n)               { return ((!std::numeric_limits<char>::is_signed) ? u.mul_unsigned_long_long(static_cast<unsigned long long>(n)) : u.mul_signed_long_long(static_cast<signed long long>(n))); }
+      inline mp_float& operator*=(mp_float& u, const signed char n)        { return u.mul_signed_long_long(n); }
+      inline mp_float& operator*=(mp_float& u, const signed short n)       { return u.mul_signed_long_long(n); }
+      inline mp_float& operator*=(mp_float& u, const signed int n)         { return u.mul_signed_long_long(n); }
+      inline mp_float& operator*=(mp_float& u, const signed long n)        { return u.mul_signed_long_long(n); }
+      inline mp_float& operator*=(mp_float& u, const signed long long n)   { return u.mul_signed_long_long(n); }
+      inline mp_float& operator*=(mp_float& u, const unsigned char n)      { return u.mul_unsigned_long_long(n); }
+      inline mp_float& operator*=(mp_float& u, const wchar_t n)            { return ((!std::numeric_limits<wchar_t>::is_signed) ? u.mul_unsigned_long_long(static_cast<unsigned long long>(n)) : u.mul_signed_long_long(static_cast<signed long long>(n))); }
+      inline mp_float& operator*=(mp_float& u, const unsigned short n)     { return u.mul_unsigned_long_long(n); }
+      inline mp_float& operator*=(mp_float& u, const unsigned int n)       { return u.mul_unsigned_long_long(n); }
+      inline mp_float& operator*=(mp_float& u, const unsigned long n)      { return u.mul_unsigned_long_long(n); }
+      inline mp_float& operator*=(mp_float& u, const unsigned long long n) { return u.mul_unsigned_long_long(n); }
+      inline mp_float& operator*=(mp_float& u, const float f)              { return u *= mp_float(f); }
+      inline mp_float& operator*=(mp_float& u, const double d)             { return u *= mp_float(d); }
+      inline mp_float& operator*=(mp_float& u, const long double ld)       { return u *= mp_float(ld); }
+
+      inline mp_float& operator/=(mp_float& u, const char n)               { return ((!std::numeric_limits<char>::is_signed) ? u.div_unsigned_long_long(static_cast<unsigned long long>(n)) : u.div_signed_long_long(static_cast<signed long long>(n))); }
+      inline mp_float& operator/=(mp_float& u, const signed char n)        { return u.div_signed_long_long(n); }
+      inline mp_float& operator/=(mp_float& u, const signed short n)       { return u.div_signed_long_long(n); }
+      inline mp_float& operator/=(mp_float& u, const signed int n)         { return u.div_signed_long_long(n); }
+      inline mp_float& operator/=(mp_float& u, const signed long n)        { return u.div_signed_long_long(n); }
+      inline mp_float& operator/=(mp_float& u, const signed long long n)   { return u.div_signed_long_long(n); }
+      inline mp_float& operator/=(mp_float& u, const unsigned char n)      { return u.div_unsigned_long_long(n); }
+      inline mp_float& operator/=(mp_float& u, const wchar_t n)            { return ((!std::numeric_limits<wchar_t>::is_signed) ? u.div_unsigned_long_long(static_cast<unsigned long long>(n)) : u.div_signed_long_long(static_cast<signed long long>(n))); }
+      inline mp_float& operator/=(mp_float& u, const unsigned short n)     { return u.div_unsigned_long_long(n); }
+      inline mp_float& operator/=(mp_float& u, const unsigned int n)       { return u.div_unsigned_long_long(n); }
+      inline mp_float& operator/=(mp_float& u, const unsigned long n)      { return u.div_unsigned_long_long(n); }
+      inline mp_float& operator/=(mp_float& u, const unsigned long long n) { return u.div_unsigned_long_long(n); }
+      inline mp_float& operator/=(mp_float& u, const float f)              { return u /= mp_float(f); }
+      inline mp_float& operator/=(mp_float& u, const double d)             { return u /= mp_float(d); }
+      inline mp_float& operator/=(mp_float& u, const long double ld)       { return u /= mp_float(ld); }
+
+      // Global comparison operators of const mp_float& with const mp_float&.
+      inline bool operator< (const mp_float& u, const mp_float& v) { return (u.cmp(v) <  static_cast<boost::int32_t>(0)); }
+      inline bool operator<=(const mp_float& u, const mp_float& v) { return (u.cmp(v) <= static_cast<boost::int32_t>(0)); }
+      inline bool operator==(const mp_float& u, const mp_float& v) { return (u.cmp(v) == static_cast<boost::int32_t>(0)); }
+      inline bool operator!=(const mp_float& u, const mp_float& v) { return (u.cmp(v) != static_cast<boost::int32_t>(0)); }
+      inline bool operator>=(const mp_float& u, const mp_float& v) { return (u.cmp(v) >= static_cast<boost::int32_t>(0)); }
+      inline bool operator> (const mp_float& u, const mp_float& v) { return (u.cmp(v) >  static_cast<boost::int32_t>(0)); }
+
+      // Global comparison operators of const mp_float& with all built-in types.
+      inline bool operator< (const mp_float& u, const char v) { return (u.cmp(mp_float(v)) <  static_cast<boost::int32_t>(0)); }
+      inline bool operator<=(const mp_float& u, const char v) { return (u.cmp(mp_float(v)) <= static_cast<boost::int32_t>(0)); }
+      inline bool operator==(const mp_float& u, const char v) { return (u.cmp(mp_float(v)) == static_cast<boost::int32_t>(0)); }
+      inline bool operator!=(const mp_float& u, const char v) { return (u.cmp(mp_float(v)) != static_cast<boost::int32_t>(0)); }
+      inline bool operator>=(const mp_float& u, const char v) { return (u.cmp(mp_float(v)) >= static_cast<boost::int32_t>(0)); }
+      inline bool operator> (const mp_float& u, const char v) { return (u.cmp(mp_float(v)) >  static_cast<boost::int32_t>(0)); }
+
+      inline bool operator< (const mp_float& u, const wchar_t v) { return (u.cmp(mp_float(v)) <  static_cast<boost::int32_t>(0)); }
+      inline bool operator<=(const mp_float& u, const wchar_t v) { return (u.cmp(mp_float(v)) <= static_cast<boost::int32_t>(0)); }
+      inline bool operator==(const mp_float& u, const wchar_t v) { return (u.cmp(mp_float(v)) == static_cast<boost::int32_t>(0)); }
+      inline bool operator!=(const mp_float& u, const wchar_t v) { return (u.cmp(mp_float(v)) != static_cast<boost::int32_t>(0)); }
+      inline bool operator>=(const mp_float& u, const wchar_t v) { return (u.cmp(mp_float(v)) >= static_cast<boost::int32_t>(0)); }
+      inline bool operator> (const mp_float& u, const wchar_t v) { return (u.cmp(mp_float(v)) >  static_cast<boost::int32_t>(0)); }
+
+      inline bool operator< (const mp_float& u, const signed char v) { return (u.cmp(mp_float(v)) <  static_cast<boost::int32_t>(0)); }
+      inline bool operator<=(const mp_float& u, const signed char v) { return (u.cmp(mp_float(v)) <= static_cast<boost::int32_t>(0)); }
+      inline bool operator==(const mp_float& u, const signed char v) { return (u.cmp(mp_float(v)) == static_cast<boost::int32_t>(0)); }
+      inline bool operator!=(const mp_float& u, const signed char v) { return (u.cmp(mp_float(v)) != static_cast<boost::int32_t>(0)); }
+      inline bool operator>=(const mp_float& u, const signed char v) { return (u.cmp(mp_float(v)) >= static_cast<boost::int32_t>(0)); }
+      inline bool operator> (const mp_float& u, const signed char v) { return (u.cmp(mp_float(v)) >  static_cast<boost::int32_t>(0)); }
+
+      inline bool operator< (const mp_float& u, const signed short v) { return (u.cmp(mp_float(v)) <  static_cast<boost::int32_t>(0)); }
+      inline bool operator<=(const mp_float& u, const signed short v) { return (u.cmp(mp_float(v)) <= static_cast<boost::int32_t>(0)); }
+      inline bool operator==(const mp_float& u, const signed short v) { return (u.cmp(mp_float(v)) == static_cast<boost::int32_t>(0)); }
+      inline bool operator!=(const mp_float& u, const signed short v) { return (u.cmp(mp_float(v)) != static_cast<boost::int32_t>(0)); }
+      inline bool operator>=(const mp_float& u, const signed short v) { return (u.cmp(mp_float(v)) >= static_cast<boost::int32_t>(0)); }
+      inline bool operator> (const mp_float& u, const signed short v) { return (u.cmp(mp_float(v)) >  static_cast<boost::int32_t>(0)); }
+
+      inline bool operator< (const mp_float& u, const signed int v) { return (u.cmp(mp_float(v)) <  static_cast<boost::int32_t>(0)); }
+      inline bool operator<=(const mp_float& u, const signed int v) { return (u.cmp(mp_float(v)) <= static_cast<boost::int32_t>(0)); }
+      inline bool operator==(const mp_float& u, const signed int v) { return (u.cmp(mp_float(v)) == static_cast<boost::int32_t>(0)); }
+      inline bool operator!=(const mp_float& u, const signed int v) { return (u.cmp(mp_float(v)) != static_cast<boost::int32_t>(0)); }
+      inline bool operator>=(const mp_float& u, const signed int v) { return (u.cmp(mp_float(v)) >= static_cast<boost::int32_t>(0)); }
+      inline bool operator> (const mp_float& u, const signed int v) { return (u.cmp(mp_float(v)) >  static_cast<boost::int32_t>(0)); }
+
+      inline bool operator< (const mp_float& u, const signed long v) { return (u.cmp(mp_float(v)) <  static_cast<boost::int32_t>(0)); }
+      inline bool operator<=(const mp_float& u, const signed long v) { return (u.cmp(mp_float(v)) <= static_cast<boost::int32_t>(0)); }
+      inline bool operator==(const mp_float& u, const signed long v) { return (u.cmp(mp_float(v)) == static_cast<boost::int32_t>(0)); }
+      inline bool operator!=(const mp_float& u, const signed long v) { return (u.cmp(mp_float(v)) != static_cast<boost::int32_t>(0)); }
+      inline bool operator>=(const mp_float& u, const signed long v) { return (u.cmp(mp_float(v)) >= static_cast<boost::int32_t>(0)); }
+      inline bool operator> (const mp_float& u, const signed long v) { return (u.cmp(mp_float(v)) >  static_cast<boost::int32_t>(0)); }
+
+      inline bool operator< (const mp_float& u, const signed long long v) { return (u.cmp(mp_float(v)) <  static_cast<boost::int32_t>(0)); }
+      inline bool operator<=(const mp_float& u, const signed long long v) { return (u.cmp(mp_float(v)) <= static_cast<boost::int32_t>(0)); }
+      inline bool operator==(const mp_float& u, const signed long long v) { return (u.cmp(mp_float(v)) == static_cast<boost::int32_t>(0)); }
+      inline bool operator!=(const mp_float& u, const signed long long v) { return (u.cmp(mp_float(v)) != static_cast<boost::int32_t>(0)); }
+      inline bool operator>=(const mp_float& u, const signed long long v) { return (u.cmp(mp_float(v)) >= static_cast<boost::int32_t>(0)); }
+      inline bool operator> (const mp_float& u, const signed long long v) { return (u.cmp(mp_float(v)) >  static_cast<boost::int32_t>(0)); }
+
+      inline bool operator< (const mp_float& u, const unsigned char v) { return (u.cmp(mp_float(v)) <  static_cast<boost::int32_t>(0)); }
+      inline bool operator<=(const mp_float& u, const unsigned char v) { return (u.cmp(mp_float(v)) <= static_cast<boost::int32_t>(0)); }
+      inline bool operator==(const mp_float& u, const unsigned char v) { return (u.cmp(mp_float(v)) == static_cast<boost::int32_t>(0)); }
+      inline bool operator!=(const mp_float& u, const unsigned char v) { return (u.cmp(mp_float(v)) != static_cast<boost::int32_t>(0)); }
+      inline bool operator>=(const mp_float& u, const unsigned char v) { return (u.cmp(mp_float(v)) >= static_cast<boost::int32_t>(0)); }
+      inline bool operator> (const mp_float& u, const unsigned char v) { return (u.cmp(mp_float(v)) >  static_cast<boost::int32_t>(0)); }
+
+      inline bool operator< (const mp_float& u, const unsigned short v) { return (u.cmp(mp_float(v)) <  static_cast<boost::int32_t>(0)); }
+      inline bool operator<=(const mp_float& u, const unsigned short v) { return (u.cmp(mp_float(v)) <= static_cast<boost::int32_t>(0)); }
+      inline bool operator==(const mp_float& u, const unsigned short v) { return (u.cmp(mp_float(v)) == static_cast<boost::int32_t>(0)); }
+      inline bool operator!=(const mp_float& u, const unsigned short v) { return (u.cmp(mp_float(v)) != static_cast<boost::int32_t>(0)); }
+      inline bool operator>=(const mp_float& u, const unsigned short v) { return (u.cmp(mp_float(v)) >= static_cast<boost::int32_t>(0)); }
+      inline bool operator> (const mp_float& u, const unsigned short v) { return (u.cmp(mp_float(v)) >  static_cast<boost::int32_t>(0)); }
+
+      inline bool operator< (const mp_float& u, const unsigned int v) { return (u.cmp(mp_float(v)) <  static_cast<boost::int32_t>(0)); }
+      inline bool operator<=(const mp_float& u, const unsigned int v) { return (u.cmp(mp_float(v)) <= static_cast<boost::int32_t>(0)); }
+      inline bool operator==(const mp_float& u, const unsigned int v) { return (u.cmp(mp_float(v)) == static_cast<boost::int32_t>(0)); }
+      inline bool operator!=(const mp_float& u, const unsigned int v) { return (u.cmp(mp_float(v)) != static_cast<boost::int32_t>(0)); }
+      inline bool operator>=(const mp_float& u, const unsigned int v) { return (u.cmp(mp_float(v)) >= static_cast<boost::int32_t>(0)); }
+      inline bool operator> (const mp_float& u, const unsigned int v) { return (u.cmp(mp_float(v)) >  static_cast<boost::int32_t>(0)); }
+
+      inline bool operator< (const mp_float& u, const unsigned long v) { return (u.cmp(mp_float(v)) <  static_cast<boost::int32_t>(0)); }
+      inline bool operator<=(const mp_float& u, const unsigned long v) { return (u.cmp(mp_float(v)) <= static_cast<boost::int32_t>(0)); }
+      inline bool operator==(const mp_float& u, const unsigned long v) { return (u.cmp(mp_float(v)) == static_cast<boost::int32_t>(0)); }
+      inline bool operator!=(const mp_float& u, const unsigned long v) { return (u.cmp(mp_float(v)) != static_cast<boost::int32_t>(0)); }
+      inline bool operator>=(const mp_float& u, const unsigned long v) { return (u.cmp(mp_float(v)) >= static_cast<boost::int32_t>(0)); }
+      inline bool operator> (const mp_float& u, const unsigned long v) { return (u.cmp(mp_float(v)) >  static_cast<boost::int32_t>(0)); }
+
+      inline bool operator< (const mp_float& u, const unsigned long long v) { return (u.cmp(mp_float(v)) <  static_cast<boost::int32_t>(0)); }
+      inline bool operator<=(const mp_float& u, const unsigned long long v) { return (u.cmp(mp_float(v)) <= static_cast<boost::int32_t>(0)); }
+      inline bool operator==(const mp_float& u, const unsigned long long v) { return (u.cmp(mp_float(v)) == static_cast<boost::int32_t>(0)); }
+      inline bool operator!=(const mp_float& u, const unsigned long long v) { return (u.cmp(mp_float(v)) != static_cast<boost::int32_t>(0)); }
+      inline bool operator>=(const mp_float& u, const unsigned long long v) { return (u.cmp(mp_float(v)) >= static_cast<boost::int32_t>(0)); }
+      inline bool operator> (const mp_float& u, const unsigned long long v) { return (u.cmp(mp_float(v)) >  static_cast<boost::int32_t>(0)); }
+
+      inline bool operator< (const mp_float& u, const float v) { return (u.cmp(mp_float(v)) <  static_cast<boost::int32_t>(0)); }
+      inline bool operator<=(const mp_float& u, const float v) { return (u.cmp(mp_float(v)) <= static_cast<boost::int32_t>(0)); }
+      inline bool operator==(const mp_float& u, const float v) { return (u.cmp(mp_float(v)) == static_cast<boost::int32_t>(0)); }
+      inline bool operator!=(const mp_float& u, const float v) { return (u.cmp(mp_float(v)) != static_cast<boost::int32_t>(0)); }
+      inline bool operator>=(const mp_float& u, const float v) { return (u.cmp(mp_float(v)) >= static_cast<boost::int32_t>(0)); }
+      inline bool operator> (const mp_float& u, const float v) { return (u.cmp(mp_float(v)) >  static_cast<boost::int32_t>(0)); }
+
+      inline bool operator< (const mp_float& u, const double v) { return (u.cmp(mp_float(v)) <  static_cast<boost::int32_t>(0)); }
+      inline bool operator<=(const mp_float& u, const double v) { return (u.cmp(mp_float(v)) <= static_cast<boost::int32_t>(0)); }
+      inline bool operator==(const mp_float& u, const double v) { return (u.cmp(mp_float(v)) == static_cast<boost::int32_t>(0)); }
+      inline bool operator!=(const mp_float& u, const double v) { return (u.cmp(mp_float(v)) != static_cast<boost::int32_t>(0)); }
+      inline bool operator>=(const mp_float& u, const double v) { return (u.cmp(mp_float(v)) >= static_cast<boost::int32_t>(0)); }
+      inline bool operator> (const mp_float& u, const double v) { return (u.cmp(mp_float(v)) >  static_cast<boost::int32_t>(0)); }
+
+      inline bool operator< (const mp_float& u, const long double v) { return (u.cmp(mp_float(v)) <  static_cast<boost::int32_t>(0)); }
+      inline bool operator<=(const mp_float& u, const long double v) { return (u.cmp(mp_float(v)) <= static_cast<boost::int32_t>(0)); }
+      inline bool operator==(const mp_float& u, const long double v) { return (u.cmp(mp_float(v)) == static_cast<boost::int32_t>(0)); }
+      inline bool operator!=(const mp_float& u, const long double v) { return (u.cmp(mp_float(v)) != static_cast<boost::int32_t>(0)); }
+      inline bool operator>=(const mp_float& u, const long double v) { return (u.cmp(mp_float(v)) >= static_cast<boost::int32_t>(0)); }
+      inline bool operator> (const mp_float& u, const long double v) { return (u.cmp(mp_float(v)) >  static_cast<boost::int32_t>(0)); }
+
+      // Global comparison operators of all built-in types with const mp_float&.
+      inline bool operator< (const char u, const mp_float& v) { return (mp_float(u).cmp(v) <  static_cast<boost::int32_t>(0)); }
+      inline bool operator<=(const char u, const mp_float& v) { return (mp_float(u).cmp(v) <= static_cast<boost::int32_t>(0)); }
+      inline bool operator==(const char u, const mp_float& v) { return (mp_float(u).cmp(v) == static_cast<boost::int32_t>(0)); }
+      inline bool operator!=(const char u, const mp_float& v) { return (mp_float(u).cmp(v) != static_cast<boost::int32_t>(0)); }
+      inline bool operator>=(const char u, const mp_float& v) { return (mp_float(u).cmp(v) >= static_cast<boost::int32_t>(0)); }
+      inline bool operator> (const char u, const mp_float& v) { return (mp_float(u).cmp(v) >  static_cast<boost::int32_t>(0)); }
+
+      inline bool operator< (const wchar_t u, const mp_float& v) { return (mp_float(u).cmp(v) <  static_cast<boost::int32_t>(0)); }
+      inline bool operator<=(const wchar_t u, const mp_float& v) { return (mp_float(u).cmp(v) <= static_cast<boost::int32_t>(0)); }
+      inline bool operator==(const wchar_t u, const mp_float& v) { return (mp_float(u).cmp(v) == static_cast<boost::int32_t>(0)); }
+      inline bool operator!=(const wchar_t u, const mp_float& v) { return (mp_float(u).cmp(v) != static_cast<boost::int32_t>(0)); }
+      inline bool operator>=(const wchar_t u, const mp_float& v) { return (mp_float(u).cmp(v) >= static_cast<boost::int32_t>(0)); }
+      inline bool operator> (const wchar_t u, const mp_float& v) { return (mp_float(u).cmp(v) >  static_cast<boost::int32_t>(0)); }
+
+      inline bool operator< (const signed char u, const mp_float& v) { return (mp_float(u).cmp(v) <  static_cast<boost::int32_t>(0)); }
+      inline bool operator<=(const signed char u, const mp_float& v) { return (mp_float(u).cmp(v) <= static_cast<boost::int32_t>(0)); }
+      inline bool operator==(const signed char u, const mp_float& v) { return (mp_float(u).cmp(v) == static_cast<boost::int32_t>(0)); }
+      inline bool operator!=(const signed char u, const mp_float& v) { return (mp_float(u).cmp(v) != static_cast<boost::int32_t>(0)); }
+      inline bool operator>=(const signed char u, const mp_float& v) { return (mp_float(u).cmp(v) >= static_cast<boost::int32_t>(0)); }
+      inline bool operator> (const signed char u, const mp_float& v) { return (mp_float(u).cmp(v) >  static_cast<boost::int32_t>(0)); }
+
+      inline bool operator< (const signed short u, const mp_float& v) { return (mp_float(u).cmp(v) <  static_cast<boost::int32_t>(0)); }
+      inline bool operator<=(const signed short u, const mp_float& v) { return (mp_float(u).cmp(v) <= static_cast<boost::int32_t>(0)); }
+      inline bool operator==(const signed short u, const mp_float& v) { return (mp_float(u).cmp(v) == static_cast<boost::int32_t>(0)); }
+      inline bool operator!=(const signed short u, const mp_float& v) { return (mp_float(u).cmp(v) != static_cast<boost::int32_t>(0)); }
+      inline bool operator>=(const signed short u, const mp_float& v) { return (mp_float(u).cmp(v) >= static_cast<boost::int32_t>(0)); }
+      inline bool operator> (const signed short u, const mp_float& v) { return (mp_float(u).cmp(v) >  static_cast<boost::int32_t>(0)); }
+
+      inline bool operator< (const signed int u, const mp_float& v) { return (mp_float(u).cmp(v) <  static_cast<boost::int32_t>(0)); }
+      inline bool operator<=(const signed int u, const mp_float& v) { return (mp_float(u).cmp(v) <= static_cast<boost::int32_t>(0)); }
+      inline bool operator==(const signed int u, const mp_float& v) { return (mp_float(u).cmp(v) == static_cast<boost::int32_t>(0)); }
+      inline bool operator!=(const signed int u, const mp_float& v) { return (mp_float(u).cmp(v) != static_cast<boost::int32_t>(0)); }
+      inline bool operator>=(const signed int u, const mp_float& v) { return (mp_float(u).cmp(v) >= static_cast<boost::int32_t>(0)); }
+      inline bool operator> (const signed int u, const mp_float& v) { return (mp_float(u).cmp(v) >  static_cast<boost::int32_t>(0)); }
+
+      inline bool operator< (const signed long u, const mp_float& v) { return (mp_float(u).cmp(v) <  static_cast<boost::int32_t>(0)); }
+      inline bool operator<=(const signed long u, const mp_float& v) { return (mp_float(u).cmp(v) <= static_cast<boost::int32_t>(0)); }
+      inline bool operator==(const signed long u, const mp_float& v) { return (mp_float(u).cmp(v) == static_cast<boost::int32_t>(0)); }
+      inline bool operator!=(const signed long u, const mp_float& v) { return (mp_float(u).cmp(v) != static_cast<boost::int32_t>(0)); }
+      inline bool operator>=(const signed long u, const mp_float& v) { return (mp_float(u).cmp(v) >= static_cast<boost::int32_t>(0)); }
+      inline bool operator> (const signed long u, const mp_float& v) { return (mp_float(u).cmp(v) >  static_cast<boost::int32_t>(0)); }
+
+      inline bool operator< (const signed long long u, const mp_float& v) { return (mp_float(u).cmp(v) <  static_cast<boost::int32_t>(0)); }
+      inline bool operator<=(const signed long long u, const mp_float& v) { return (mp_float(u).cmp(v) <= static_cast<boost::int32_t>(0)); }
+      inline bool operator==(const signed long long u, const mp_float& v) { return (mp_float(u).cmp(v) == static_cast<boost::int32_t>(0)); }
+      inline bool operator!=(const signed long long u, const mp_float& v) { return (mp_float(u).cmp(v) != static_cast<boost::int32_t>(0)); }
+      inline bool operator>=(const signed long long u, const mp_float& v) { return (mp_float(u).cmp(v) >= static_cast<boost::int32_t>(0)); }
+      inline bool operator> (const signed long long u, const mp_float& v) { return (mp_float(u).cmp(v) >  static_cast<boost::int32_t>(0)); }
+
+      inline bool operator< (const unsigned char u, const mp_float& v) { return (mp_float(u).cmp(v) <  static_cast<boost::int32_t>(0)); }
+      inline bool operator<=(const unsigned char u, const mp_float& v) { return (mp_float(u).cmp(v) <= static_cast<boost::int32_t>(0)); }
+      inline bool operator==(const unsigned char u, const mp_float& v) { return (mp_float(u).cmp(v) == static_cast<boost::int32_t>(0)); }
+      inline bool operator!=(const unsigned char u, const mp_float& v) { return (mp_float(u).cmp(v) != static_cast<boost::int32_t>(0)); }
+      inline bool operator>=(const unsigned char u, const mp_float& v) { return (mp_float(u).cmp(v) >= static_cast<boost::int32_t>(0)); }
+      inline bool operator> (const unsigned char u, const mp_float& v) { return (mp_float(u).cmp(v) >  static_cast<boost::int32_t>(0)); }
+
+      inline bool operator< (const unsigned short u, const mp_float& v) { return (mp_float(u).cmp(v) <  static_cast<boost::int32_t>(0)); }
+      inline bool operator<=(const unsigned short u, const mp_float& v) { return (mp_float(u).cmp(v) <= static_cast<boost::int32_t>(0)); }
+      inline bool operator==(const unsigned short u, const mp_float& v) { return (mp_float(u).cmp(v) == static_cast<boost::int32_t>(0)); }
+      inline bool operator!=(const unsigned short u, const mp_float& v) { return (mp_float(u).cmp(v) != static_cast<boost::int32_t>(0)); }
+      inline bool operator>=(const unsigned short u, const mp_float& v) { return (mp_float(u).cmp(v) >= static_cast<boost::int32_t>(0)); }
+      inline bool operator> (const unsigned short u, const mp_float& v) { return (mp_float(u).cmp(v) >  static_cast<boost::int32_t>(0)); }
+
+      inline bool operator< (const unsigned int u, const mp_float& v) { return (mp_float(u).cmp(v) <  static_cast<boost::int32_t>(0)); }
+      inline bool operator<=(const unsigned int u, const mp_float& v) { return (mp_float(u).cmp(v) <= static_cast<boost::int32_t>(0)); }
+      inline bool operator==(const unsigned int u, const mp_float& v) { return (mp_float(u).cmp(v) == static_cast<boost::int32_t>(0)); }
+      inline bool operator!=(const unsigned int u, const mp_float& v) { return (mp_float(u).cmp(v) != static_cast<boost::int32_t>(0)); }
+      inline bool operator>=(const unsigned int u, const mp_float& v) { return (mp_float(u).cmp(v) >= static_cast<boost::int32_t>(0)); }
+      inline bool operator> (const unsigned int u, const mp_float& v) { return (mp_float(u).cmp(v) >  static_cast<boost::int32_t>(0)); }
+
+      inline bool operator< (const unsigned long u, const mp_float& v) { return (mp_float(u).cmp(v) <  static_cast<boost::int32_t>(0)); }
+      inline bool operator<=(const unsigned long u, const mp_float& v) { return (mp_float(u).cmp(v) <= static_cast<boost::int32_t>(0)); }
+      inline bool operator==(const unsigned long u, const mp_float& v) { return (mp_float(u).cmp(v) == static_cast<boost::int32_t>(0)); }
+      inline bool operator!=(const unsigned long u, const mp_float& v) { return (mp_float(u).cmp(v) != static_cast<boost::int32_t>(0)); }
+      inline bool operator>=(const unsigned long u, const mp_float& v) { return (mp_float(u).cmp(v) >= static_cast<boost::int32_t>(0)); }
+      inline bool operator> (const unsigned long u, const mp_float& v) { return (mp_float(u).cmp(v) >  static_cast<boost::int32_t>(0)); }
+
+      inline bool operator< (const unsigned long long u, const mp_float& v) { return (mp_float(u).cmp(v) <  static_cast<boost::int32_t>(0)); }
+      inline bool operator<=(const unsigned long long u, const mp_float& v) { return (mp_float(u).cmp(v) <= static_cast<boost::int32_t>(0)); }
+      inline bool operator==(const unsigned long long u, const mp_float& v) { return (mp_float(u).cmp(v) == static_cast<boost::int32_t>(0)); }
+      inline bool operator!=(const unsigned long long u, const mp_float& v) { return (mp_float(u).cmp(v) != static_cast<boost::int32_t>(0)); }
+      inline bool operator>=(const unsigned long long u, const mp_float& v) { return (mp_float(u).cmp(v) >= static_cast<boost::int32_t>(0)); }
+      inline bool operator> (const unsigned long long u, const mp_float& v) { return (mp_float(u).cmp(v) >  static_cast<boost::int32_t>(0)); }
+
+      inline bool operator< (const float u, const mp_float& v) { return (mp_float(u).cmp(v) <  static_cast<boost::int32_t>(0)); }
+      inline bool operator<=(const float u, const mp_float& v) { return (mp_float(u).cmp(v) <= static_cast<boost::int32_t>(0)); }
+      inline bool operator==(const float u, const mp_float& v) { return (mp_float(u).cmp(v) == static_cast<boost::int32_t>(0)); }
+      inline bool operator!=(const float u, const mp_float& v) { return (mp_float(u).cmp(v) != static_cast<boost::int32_t>(0)); }
+      inline bool operator>=(const float u, const mp_float& v) { return (mp_float(u).cmp(v) >= static_cast<boost::int32_t>(0)); }
+      inline bool operator> (const float u, const mp_float& v) { return (mp_float(u).cmp(v) >  static_cast<boost::int32_t>(0)); }
+
+      inline bool operator< (const double u, const mp_float& v) { return (mp_float(u).cmp(v) <  static_cast<boost::int32_t>(0)); }
+      inline bool operator<=(const double u, const mp_float& v) { return (mp_float(u).cmp(v) <= static_cast<boost::int32_t>(0)); }
+      inline bool operator==(const double u, const mp_float& v) { return (mp_float(u).cmp(v) == static_cast<boost::int32_t>(0)); }
+      inline bool operator!=(const double u, const mp_float& v) { return (mp_float(u).cmp(v) != static_cast<boost::int32_t>(0)); }
+      inline bool operator>=(const double u, const mp_float& v) { return (mp_float(u).cmp(v) >= static_cast<boost::int32_t>(0)); }
+      inline bool operator> (const double u, const mp_float& v) { return (mp_float(u).cmp(v) >  static_cast<boost::int32_t>(0)); }
+
+      inline bool operator< (const long double u, const mp_float& v) { return (mp_float(u).cmp(v) <  static_cast<boost::int32_t>(0)); }
+      inline bool operator<=(const long double u, const mp_float& v) { return (mp_float(u).cmp(v) <= static_cast<boost::int32_t>(0)); }
+      inline bool operator==(const long double u, const mp_float& v) { return (mp_float(u).cmp(v) == static_cast<boost::int32_t>(0)); }
+      inline bool operator!=(const long double u, const mp_float& v) { return (mp_float(u).cmp(v) != static_cast<boost::int32_t>(0)); }
+      inline bool operator>=(const long double u, const mp_float& v) { return (mp_float(u).cmp(v) >= static_cast<boost::int32_t>(0)); }
+      inline bool operator> (const long double u, const mp_float& v) { return (mp_float(u).cmp(v) >  static_cast<boost::int32_t>(0)); }
+    }
+  }
+
+#endif // _MP_FLOAT_2004_02_08_HPP_

Added: sandbox/multiprecision/boost/mp_float_base.hpp
==============================================================================
--- (empty file)
+++ sandbox/multiprecision/boost/mp_float_base.hpp	2011-09-21 17:00:15 EDT (Wed, 21 Sep 2011)
@@ -0,0 +1,277 @@
+
+//          Copyright Christopher Kormanyos 2002 - 2011.
+// Distributed under the Boost Software License, Version 1.0.
+//    (See accompanying file LICENSE_1_0.txt or copy at
+//          http://www.boost.org/LICENSE_1_0.txt)
+
+// This work is based on an earlier work:
+// "Algorithm 910: A Portable C++ Multiple-Precision System for Special-Function Calculations",
+// in ACM TOMS, {VOL 37, ISSUE 4, (February 2011)} (C) ACM, 2011. http://doi.acm.org/10.1145/1916461.1916469
+
+#ifndef _MP_FLOAT_BASE_2004_02_08_HPP_
+  #define _MP_FLOAT_BASE_2004_02_08_HPP_
+
+  #include <iostream>
+  #include <limits>
+  #include <string>
+  #include <cmath>
+
+  #include <boost/cstdint.hpp>
+
+  // Select the number of decimal digits in mp_float
+  // by setting the value of BOOST_MULTIPRECISION_BACKEND_MP_FLOAT_DIGITS10.
+  // The supported range is 30-300.
+  // Note: This is a compile-time constant.
+
+  #if !defined(BOOST_MULTIPRECISION_BACKEND_MP_FLOAT_DIGITS10)
+    #define BOOST_MULTIPRECISION_BACKEND_MP_FLOAT_DIGITS10 50
+  #endif
+  #define BOOST_MULTIPRECISION_BACKEND_MP_FLOAT_DIGITS10_LIMIT 300
+
+  namespace boost
+  {
+    namespace multiprecision
+    {
+    #if defined(BOOST_MULTIPRECISION_BACKEND_MP_FLOAT_TYPE_EFX)
+      class mp_float_efx;
+      typedef mp_float_efx mp_float;
+    #elif defined(BOOST_MULTIPRECISION_BACKEND_MP_FLOAT_TYPE_GMP)
+      class mp_float_gmp;
+      typedef mp_float_gmp mp_float;
+    #elif defined(BOOST_MULTIPRECISION_BACKEND_MP_FLOAT_TYPE_MPFR)
+      class mp_float_mpfr;
+      typedef mp_float_mpfr mp_float;
+    #else
+      #error The MP float type is undefined! Define the mp_float type!
+    #endif
+    }
+  }
+
+  namespace boost
+  {
+    namespace multiprecision
+    {
+      class mp_float_base
+      {
+      public:
+        static const boost::int32_t mp_digits10_setting = BOOST_MULTIPRECISION_BACKEND_MP_FLOAT_DIGITS10;
+        static const boost::int32_t mp_digits10_limit   = BOOST_MULTIPRECISION_BACKEND_MP_FLOAT_DIGITS10_LIMIT;
+        static const boost::int32_t mp_digits10         = ((mp_digits10_setting < static_cast<boost::int32_t>(30)) ? static_cast<boost::int32_t>(30) : ((mp_digits10_setting > mp_digits10_limit) ? mp_digits10_limit : mp_digits10_setting));
+        static const boost::int32_t mp_digits10_extra   = static_cast<boost::int32_t>(((static_cast<boost::int64_t>(mp_digits10) * 15LL) + 50LL) / 100LL);
+        static const boost::int32_t mp_max_digits10     = static_cast<boost::int32_t>(mp_digits10 + ((mp_digits10_extra < static_cast<boost::int32_t>(5)) ? static_cast<boost::int32_t>(5) : ((mp_digits10_extra > static_cast<boost::int32_t>(30)) ? static_cast<boost::int32_t>(30) : mp_digits10_extra)));
+
+        virtual ~mp_float_base() { }
+
+        // Specific special values.
+        virtual const mp_float& my_value_nan(void) const = 0;
+        virtual const mp_float& my_value_inf(void) const = 0;
+        const mp_float& my_value_max(void) const;
+        const mp_float& my_value_min(void) const;
+
+        virtual void precision(const boost::int32_t) = 0;
+
+        // Basic operations.
+        virtual mp_float& operator= (const mp_float&) = 0;
+        virtual mp_float& operator+=(const mp_float&) = 0;
+        virtual mp_float& operator-=(const mp_float&) = 0;
+        virtual mp_float& operator*=(const mp_float&) = 0;
+        virtual mp_float& operator/=(const mp_float&) = 0;
+        virtual mp_float& add_unsigned_long_long(const unsigned long long) = 0;
+        virtual mp_float& sub_unsigned_long_long(const unsigned long long) = 0;
+        virtual mp_float& mul_unsigned_long_long(const unsigned long long) = 0;
+        virtual mp_float& div_unsigned_long_long(const unsigned long long) = 0;
+
+        mp_float& add_signed_long_long(const signed long long);
+        mp_float& sub_signed_long_long(const signed long long);
+        mp_float& mul_signed_long_long(const signed long long);
+        mp_float& div_signed_long_long(const signed long long);
+
+        // Elementary primitives.
+        virtual mp_float& calculate_inv (void) = 0;
+        virtual mp_float& calculate_sqrt(void) = 0;
+        virtual mp_float& negate(void) = 0;
+
+        // Comparison functions.
+        virtual boost::int32_t cmp(const mp_float&) const = 0;
+
+        virtual bool isnan   (void) const = 0;
+        virtual bool isinf   (void) const = 0;
+        virtual bool isfinite(void) const = 0;
+
+        virtual bool iszero(void) const = 0;
+        virtual bool isone (void) const = 0;
+        virtual bool isint (void) const = 0;
+        virtual bool isneg (void) const = 0;
+                bool ispos (void) const { return (!isneg()); }
+
+        // Operators pre-increment and pre-decrement.
+        virtual mp_float_base& operator++(void) = 0;
+        virtual mp_float_base& operator--(void) = 0;
+
+        // Fast order-10 range check.
+        boost::int64_t order(void) const { return get_order_fast(); }
+
+        // Conversion routines.
+        virtual void               extract_parts             (double&, boost::int64_t&) const = 0;
+        virtual double             extract_double            (void) const = 0;
+        virtual long double        extract_long_double       (void) const = 0;
+        virtual signed long long   extract_signed_long_long  (void) const = 0;
+        virtual unsigned long long extract_unsigned_long_long(void) const = 0;
+        virtual mp_float            extract_integer_part     (void) const = 0;
+        virtual mp_float            extract_decimal_part     (void) const = 0;
+
+        // Formated Output routine.
+        void wr_string(std::string& str, std::ostream& os) const;
+        virtual bool rd_string(const char* const) = 0;
+
+        // Cast operators with all built-in types.
+        operator char() const               { return (std::numeric_limits<char>::is_signed ? static_cast<char>   (extract_signed_long_long()) : static_cast<char>   (extract_unsigned_long_long())); }
+        operator wchar_t() const            { return (std::numeric_limits<char>::is_signed ? static_cast<wchar_t>(extract_signed_long_long()) : static_cast<wchar_t>(extract_unsigned_long_long())); }
+        operator signed char() const        { return static_cast<signed char>       (extract_signed_long_long()); }
+        operator signed short() const       { return static_cast<signed short>      (extract_signed_long_long()); }
+        operator signed int() const         { return static_cast<signed int>        (extract_signed_long_long()); }
+        operator signed long() const        { return static_cast<signed long>       (extract_signed_long_long()); }
+        operator signed long long() const   { return static_cast<signed long long>  (extract_signed_long_long()); }
+        operator unsigned char() const      { return static_cast<unsigned char>     (extract_unsigned_long_long()); }
+        operator unsigned short() const     { return static_cast<unsigned short>    (extract_unsigned_long_long()); }
+        operator unsigned int() const       { return static_cast<unsigned int>      (extract_unsigned_long_long()); }
+        operator unsigned long() const      { return static_cast<unsigned long>     (extract_unsigned_long_long()); }
+        operator unsigned long long() const { return static_cast<unsigned long long>(extract_unsigned_long_long()); }
+        operator float() const              { return static_cast<float>(extract_double()); }
+        operator double() const             { return extract_double(); }
+        operator long double() const        { return extract_long_double(); }
+
+        // Specific higher functions which might be present in the MP implementation.
+        virtual bool has_its_own_ldexp        (void) const { return false; }
+        virtual bool has_its_own_frexp        (void) const { return false; }
+        virtual bool has_its_own_fmod         (void) const { return false; }
+        virtual bool has_its_own_cbrt         (void) const { return false; }
+        virtual bool has_its_own_rootn        (void) const { return false; }
+        virtual bool has_its_own_exp          (void) const { return false; }
+        virtual bool has_its_own_log          (void) const { return false; }
+        virtual bool has_its_own_sin          (void) const { return false; }
+        virtual bool has_its_own_cos          (void) const { return false; }
+        virtual bool has_its_own_tan          (void) const { return false; }
+        virtual bool has_its_own_asin         (void) const { return false; }
+        virtual bool has_its_own_acos         (void) const { return false; }
+        virtual bool has_its_own_atan         (void) const { return false; }
+        virtual bool has_its_own_sinh         (void) const { return false; }
+        virtual bool has_its_own_cosh         (void) const { return false; }
+        virtual bool has_its_own_tanh         (void) const { return false; }
+        virtual bool has_its_own_asinh        (void) const { return false; }
+        virtual bool has_its_own_acosh        (void) const { return false; }
+        virtual bool has_its_own_atanh        (void) const { return false; }
+        virtual bool has_its_own_gamma        (void) const { return false; }
+        virtual bool has_its_own_riemann_zeta (void) const { return false; }
+        virtual bool has_its_own_cyl_bessel_jn(void) const { return false; }
+        virtual bool has_its_own_cyl_bessel_yn(void) const { return false; }
+
+        static mp_float my_ldexp        (const mp_float&, int);
+        static mp_float my_frexp        (const mp_float&, int*);
+        static mp_float my_fmod         (const mp_float&, const mp_float&);
+        static mp_float my_cbrt         (const mp_float&);
+        static mp_float my_rootn        (const mp_float&, const boost::uint32_t);
+        static mp_float my_exp          (const mp_float&);
+        static mp_float my_log          (const mp_float&);
+        static mp_float my_sin          (const mp_float&);
+        static mp_float my_cos          (const mp_float&);
+        static mp_float my_tan          (const mp_float&);
+        static mp_float my_asin         (const mp_float&);
+        static mp_float my_acos         (const mp_float&);
+        static mp_float my_atan         (const mp_float&);
+        static mp_float my_sinh         (const mp_float&);
+        static mp_float my_cosh         (const mp_float&);
+        static mp_float my_tanh         (const mp_float&);
+        static mp_float my_asinh        (const mp_float&);
+        static mp_float my_acosh        (const mp_float&);
+        static mp_float my_atanh        (const mp_float&);
+        static mp_float my_gamma        (const mp_float&);
+        static mp_float my_riemann_zeta (const mp_float&);
+        static mp_float my_cyl_bessel_jn(const boost::int32_t, const mp_float&);
+        static mp_float my_cyl_bessel_yn(const boost::int32_t, const mp_float&);
+
+        static bool char_is_nonzero_predicate(const char& c) { return (c != static_cast<char>('0')); }
+
+      protected:
+        inline mp_float_base();
+
+        // Emphasize: This template class can be used with native floating-point
+        // types like float, double and long double. Note: For long double,
+        // you need to verify that the mantissa fits in unsigned long long.
+        template<typename native_float_type>
+        class native_float_parts
+        {
+        public:
+          native_float_parts(const native_float_type f) : u(0uLL), e(0) { make_parts(f); }
+
+          const unsigned long long& get_mantissa(void) const { return u; }
+          const int& get_exponent(void) const { return e; }
+
+        private:
+          native_float_parts();
+          native_float_parts(const native_float_parts&);
+          const native_float_parts& operator=(const native_float_parts&);
+
+          unsigned long long u;
+          int e;
+
+          void make_parts(const native_float_type f)
+          {
+            if(f == static_cast<native_float_type>(0.0)) { return; }
+
+            // Get the fraction and base-2 exponent.
+            native_float_type man = ::frexp(f, &e);
+
+            boost::uint32_t n2 = 0u;
+
+            for(boost::uint32_t i = static_cast<boost::uint32_t>(0u); i < static_cast<boost::uint32_t>(std::numeric_limits<native_float_type>::digits); i++)
+            {
+              // Extract the mantissa of the floating-point type in base-2
+              // (yes, one bit at a time) and store it in an unsigned long long.
+              // TBD: Is this really portable?
+              man *= 2;
+
+              n2   = static_cast<boost::uint32_t>(man);
+              man -= static_cast<native_float_type>(n2);
+
+              if(n2 != static_cast<boost::uint32_t>(0u))
+              {
+                u |= 1u;
+              }
+
+              if(i < static_cast<boost::uint32_t>(std::numeric_limits<native_float_type>::digits - 1))
+              {
+                u <<= 1u;
+              }
+            }
+
+            // Ensure that the value is normalized and adjust the exponent.
+            u |= static_cast<unsigned long long>(1uLL << (std::numeric_limits<native_float_type>::digits - 1));
+            e -= 1;
+          }
+        };
+
+      private:
+        static bool digits_match_lib_dll_is_ok;
+
+        virtual boost::int64_t get_order_exact(void) const = 0;
+        virtual boost::int64_t get_order_fast(void) const = 0;
+        virtual void get_output_string(std::string& str, boost::int64_t& my_exp, const std::size_t number_of_digits) const = 0;
+      };
+
+      // Create a loud link error if the digits in the
+      // mp_float headers mismatch those in a Lib or DLL.
+      template<const boost::int32_t digits10> boost::int32_t digits_match_lib_dll(void);        // There's no function body here.
+      template<> boost::int32_t digits_match_lib_dll<mp_float_base::mp_digits10>(void); // The function body is in mp_float_base.cpp.
+
+      inline mp_float_base::mp_float_base()
+      {
+        digits_match_lib_dll_is_ok = (digits_match_lib_dll<mp_digits10>() == mp_digits10);
+      }
+
+      std::ostream& operator<<(std::ostream& os, const mp_float_base& f);
+      std::istream& operator>>(std::istream& is, mp_float_base& f);
+    }
+  }
+
+#endif // _MP_FLOAT_BASE_2004_02_08_HPP_

Added: sandbox/multiprecision/boost/mp_float_efx.hpp
==============================================================================
--- (empty file)
+++ sandbox/multiprecision/boost/mp_float_efx.hpp	2011-09-21 17:00:15 EDT (Wed, 21 Sep 2011)
@@ -0,0 +1,170 @@
+
+//          Copyright Christopher Kormanyos 2002 - 2011.
+// Distributed under the Boost Software License, Version 1.0.
+//    (See accompanying file LICENSE_1_0.txt or copy at
+//          http://www.boost.org/LICENSE_1_0.txt)
+
+// This work is based on an earlier work:
+// "Algorithm 910: A Portable C++ Multiple-Precision System for Special-Function Calculations",
+// in ACM TOMS, {VOL 37, ISSUE 4, (February 2011)} (C) ACM, 2011. http://doi.acm.org/10.1145/1916461.1916469
+
+#ifndef _MP_FLOAT_EFX_2004_02_08_HPP_
+  #define _MP_FLOAT_EFX_2004_02_08_HPP_
+
+  #if defined(__INTEL_COMPILER)
+    #pragma warning (disable:981)
+  #endif
+
+  #include <cmath>
+  #include <string>
+
+  #include <boost/array.hpp>
+  #include <boost/multiprecision/mp_float_base.hpp>
+
+  namespace boost
+  {
+    namespace multiprecision
+    {
+      class mp_float_efx : public mp_float_base
+      {
+      public:
+        static const boost::int32_t mp_radix  = static_cast<boost::int32_t>(10);
+        static const boost::int32_t mp_digits = mp_digits10;
+
+        static const boost::int64_t mp_max_exp   = static_cast<boost::int64_t>(+9223372036854775795LL);
+        static const boost::int64_t mp_min_exp   = static_cast<boost::int64_t>(-9223372036854775795LL);
+        static const boost::int64_t mp_max_exp10 = static_cast<boost::int64_t>(+3063937869882635616LL); // Approx. [mp_max_exp / log10(2)], also an even multiple of 8
+        static const boost::int64_t mp_min_exp10 = static_cast<boost::int64_t>(-3063937869882635616LL);
+
+        static const boost::int32_t mp_elem_digits10 = static_cast<boost::int32_t>(8);
+
+      private:
+        static const boost::int32_t mp_digits10_num_base = static_cast<boost::int32_t>((mp_max_digits10 / mp_elem_digits10) + (((mp_max_digits10 % mp_elem_digits10) != 0) ? 1 : 0));
+        static const boost::int32_t mp_elem_number       = static_cast<boost::int32_t>(mp_digits10_num_base + 2);
+
+        typedef enum enum_fpclass
+        {
+          mp_finite,
+          mp_inf,
+          mp_NaN
+        }
+        t_fpclass;
+
+        static const boost::int32_t mp_elem_mask = static_cast<boost::int32_t>(100000000);
+
+        typedef boost::array<boost::uint32_t, static_cast<std::size_t>(mp_elem_number)> array_type;
+
+        array_type data;
+        boost::int64_t      exp;
+        bool       neg;
+        t_fpclass  fpclass;
+        boost::int32_t      prec_elem;
+
+      public:
+        // Constructors
+        mp_float_efx() : data     (),
+                         exp      (static_cast<boost::int64_t>(0)),
+                         neg      (false),
+                         fpclass  (mp_finite),
+                         prec_elem(mp_elem_number) { }
+
+        mp_float_efx(const char n);
+        mp_float_efx(const signed char n);
+        mp_float_efx(const unsigned char n);
+        mp_float_efx(const wchar_t n);
+        mp_float_efx(const signed short n);
+        mp_float_efx(const unsigned short n);
+        mp_float_efx(const signed int n);
+        mp_float_efx(const unsigned int n);
+        mp_float_efx(const signed long n);
+        mp_float_efx(const unsigned long n);
+        mp_float_efx(const signed long long n);
+        mp_float_efx(const unsigned long long n);
+        mp_float_efx(const float f);
+        mp_float_efx(const double d);
+        mp_float_efx(const long double ld);
+        mp_float_efx(const char* const s);
+        mp_float_efx(const std::string& str);
+
+        mp_float_efx(const mp_float_efx& f) : data     (f.data),
+                                              exp      (f.exp),
+                                              neg      (f.neg),
+                                              fpclass  (f.fpclass),
+                                              prec_elem(f.prec_elem) { }
+
+        // Constructor from mantissa and exponent.
+        mp_float_efx(const double mantissa, const boost::int64_t exponent);
+
+        virtual ~mp_float_efx() { }
+
+      public:
+        virtual boost::int32_t cmp(const mp_float_efx& v) const;
+
+        // Specific special values.
+        virtual const mp_float_efx& my_value_nan(void) const;
+        virtual const mp_float_efx& my_value_inf(void) const;
+
+        virtual void precision(const boost::int32_t prec_digits);
+
+        // Basic operations.
+        virtual mp_float_efx& operator= (const mp_float_efx& v);
+        virtual mp_float_efx& operator+=(const mp_float_efx& v);
+        virtual mp_float_efx& operator-=(const mp_float_efx& v);
+        virtual mp_float_efx& operator*=(const mp_float_efx& v);
+        virtual mp_float_efx& operator/=(const mp_float_efx& v);
+        virtual mp_float_efx& add_unsigned_long_long(const unsigned long long n);
+        virtual mp_float_efx& sub_unsigned_long_long(const unsigned long long n);
+        virtual mp_float_efx& mul_unsigned_long_long(const unsigned long long n);
+        virtual mp_float_efx& div_unsigned_long_long(const unsigned long long n);
+
+        // Elementary primitives.
+        virtual mp_float_efx& calculate_inv (void);
+        virtual mp_float_efx& calculate_sqrt(void);
+        virtual mp_float_efx& negate(void) { if(!iszero()) { neg = (!neg); } return *this; }
+
+        // Comparison functions
+        virtual bool isnan   (void) const { return (fpclass == mp_NaN); }
+        virtual bool isinf   (void) const { return (fpclass == mp_inf); }
+        virtual bool isfinite(void) const { return (fpclass == mp_finite); }
+
+        virtual bool iszero (void) const;
+        virtual bool isone  (void) const;
+        virtual bool isint  (void) const;
+        virtual bool isneg  (void) const { return neg; }
+
+        // Operators pre-increment and pre-decrement
+        virtual mp_float_efx& operator++(void);
+        virtual mp_float_efx& operator--(void);
+
+        // Conversion routines
+        virtual void               extract_parts             (double& mantissa, boost::int64_t& exponent) const;
+        virtual double             extract_double            (void) const;
+        virtual long double        extract_long_double       (void) const;
+        virtual signed long long   extract_signed_long_long  (void) const;
+        virtual unsigned long long extract_unsigned_long_long(void) const;
+        virtual mp_float_efx       extract_integer_part      (void) const;
+        virtual mp_float_efx       extract_decimal_part      (void) const;
+
+      private:
+        static bool data_elem_is_non_zero_predicate(const boost::uint32_t& d) { return (d != static_cast<boost::uint32_t>(0u)); }
+        static bool data_elem_is_non_nine_predicate(const boost::uint32_t& d) { return (d != static_cast<boost::uint32_t>(mp_float_efx::mp_elem_mask - 1)); }
+
+        void from_unsigned_long_long(const unsigned long long u);
+        void from_unsigned_long(const unsigned long u);
+
+        boost::int32_t cmp_data(const array_type& vd) const;
+
+        static void mul_loop_uv(const boost::uint32_t* const u, const boost::uint32_t* const v, boost::uint32_t* const w, const boost::int32_t p);
+        static boost::uint32_t mul_loop_n(boost::uint32_t* const u, boost::uint32_t n, const boost::int32_t p);
+        static boost::uint32_t div_loop_n(boost::uint32_t* const u, boost::uint32_t n, const boost::int32_t p);
+
+        virtual boost::int64_t get_order_exact(void) const { return get_order_fast(); }
+        virtual boost::int64_t get_order_fast (void) const;
+        virtual void get_output_string(std::string& str, boost::int64_t& my_exp, const std::size_t number_of_digits) const;
+
+        virtual bool rd_string(const char* const s);
+      };
+    }
+  }
+
+#endif // _MP_FLOAT_EFX_2004_02_08_HPP_

Added: sandbox/multiprecision/boost/mp_float_functions.hpp
==============================================================================
--- (empty file)
+++ sandbox/multiprecision/boost/mp_float_functions.hpp	2011-09-21 17:00:15 EDT (Wed, 21 Sep 2011)
@@ -0,0 +1,294 @@
+
+//          Copyright Christopher Kormanyos 2002 - 2011.
+// Distributed under the Boost Software License, Version 1.0.
+//    (See accompanying file LICENSE_1_0.txt or copy at
+//          http://www.boost.org/LICENSE_1_0.txt)
+
+// This work is based on an earlier work:
+// "Algorithm 910: A Portable C++ Multiple-Precision System for Special-Function Calculations",
+// in ACM TOMS, {VOL 37, ISSUE 4, (February 2011)} (C) ACM, 2011. http://doi.acm.org/10.1145/1916461.1916469
+
+#ifndef _MP_FLOAT_FUNCTIONS_HPP_
+  #define _MP_FLOAT_FUNCTIONS_HPP_
+
+  #include <vector>
+  #include <deque>
+
+  #include <boost/multiprecision/mp_float.hpp>
+  #include <boost/multiprecision/mp_complex.hpp>
+
+  namespace boost
+  {
+    namespace multiprecision
+    {
+      inline boost::int32_t   max_iteration(void)    { return static_cast<boost::int32_t>(10000); }
+      inline boost::int64_t   tol (void)             { return static_cast<boost::int64_t>(boost::multiprecision::mp_float::mp_max_digits10); }
+
+      inline boost::multiprecision::mp_float fabs(const boost::multiprecision::mp_float& x) { return (x.isneg() ? boost::multiprecision::mp_float(x).negate() : x); }
+
+      boost::multiprecision::mp_float floor(const boost::multiprecision::mp_float& x);
+      boost::multiprecision::mp_float ceil (const boost::multiprecision::mp_float& x);
+      boost::int32_t   sgn   (const boost::multiprecision::mp_float& x);
+      boost::multiprecision::mp_float ldexp(const boost::multiprecision::mp_float& v, int e);
+      boost::multiprecision::mp_float frexp(const boost::multiprecision::mp_float& v, int* expon);
+      boost::multiprecision::mp_float fmod (const boost::multiprecision::mp_float& v1, const boost::multiprecision::mp_float& v2);
+
+             bool isnan(const double x);
+      inline bool isnan(const boost::multiprecision::mp_float& x)       { return x.isnan(); }
+      inline bool isnan(const boost::multiprecision::mp_complex& z)     { return z.isnan(); }
+
+             bool isfinite(const double x);
+      inline bool isfinite(const boost::multiprecision::mp_float& x)    { return x.isfinite(); }
+      inline bool isfinite(const boost::multiprecision::mp_complex& z)  { return z.isfinite(); }
+
+             bool isinf(const double x);
+      inline bool isinf(const boost::multiprecision::mp_float& x)       { return x.isinf(); }
+      inline bool isinf(const boost::multiprecision::mp_complex& z)     { return z.isinf(); }
+
+      inline bool isneg(const double x)          { return (x < 0.0); }
+      inline bool isneg(const boost::multiprecision::mp_float& x)       { return x.isneg(); }
+      inline bool isneg(const boost::multiprecision::mp_complex& z)     { return z.isneg(); }
+
+      inline boost::multiprecision::mp_float abs (const boost::multiprecision::mp_float& x)      { return boost::multiprecision::fabs(x); }
+      inline boost::multiprecision::mp_float real(const boost::multiprecision::mp_float& x)      { return x; }
+      inline boost::multiprecision::mp_float imag(const boost::multiprecision::mp_float&)        { return boost::multiprecision::zero(); }
+
+      inline bool ispos(const double x)          { return !isneg(x); }
+      inline bool ispos(const boost::multiprecision::mp_float& x)       { return !x.isneg(); }
+      inline bool ispos(const boost::multiprecision::mp_complex& z)     { return !z.isneg(); }
+
+             bool isint(const double x);
+      inline bool isint(const boost::multiprecision::mp_float& x)       { return x.isint(); }
+      inline bool isint(const boost::multiprecision::mp_complex& z)     { return z.isint(); }
+
+      inline bool isone(const double x)                                 { return (::fabs(1.0 - x) < (std::numeric_limits<double>::min)() * 2); }
+      inline bool isone(const boost::multiprecision::mp_float& x)       { return x.isone(); }
+      inline bool isone(const boost::multiprecision::mp_complex& z)     { return z.isone(); }
+
+      inline bool iszero(const double x)                                { return (::fabs(x) < (std::numeric_limits<double>::min)() * 2); }
+      inline bool iszero(const boost::multiprecision::mp_float& x)      { return x.iszero(); }
+      inline bool iszero(const boost::multiprecision::mp_complex& z)    { return z.iszero(); }
+
+      boost::multiprecision::mp_float integer_part(const boost::multiprecision::mp_float& x);
+      boost::multiprecision::mp_float decimal_part(const boost::multiprecision::mp_float& x);
+
+      void to_parts(const boost::multiprecision::mp_float& x, double& mantissa, boost::int64_t& exponent);
+
+      inline double to_double(const double& x)   { return x; }
+             double to_double(const boost::multiprecision::mp_float& x);
+             double to_double(const boost::multiprecision::mp_complex& z);
+
+      inline boost::int64_t order_of(const double x)    { return static_cast<boost::int64_t>(static_cast<boost::int32_t>(::log10(x))); }
+      inline boost::int64_t order_of(const boost::multiprecision::mp_float& x) { return x.order(); }
+
+      boost::int64_t to_int64(const double x);
+      boost::int64_t to_int64(const boost::multiprecision::mp_float& x);
+      boost::int64_t to_int64(const boost::multiprecision::mp_complex& z);
+
+      boost::int32_t to_int32(const double x);
+      boost::int32_t to_int32(const boost::multiprecision::mp_float& x);
+      boost::int32_t to_int32(const boost::multiprecision::mp_complex& z);
+
+      bool small_arg(const double x);
+      bool small_arg(const boost::multiprecision::mp_float& x);
+      bool small_arg(const boost::multiprecision::mp_complex& z);
+
+      bool large_arg(const double x);
+      bool large_arg(const boost::multiprecision::mp_float& x);
+      bool large_arg(const boost::multiprecision::mp_complex& z);
+
+      bool near_one(const double x);
+      bool near_one(const boost::multiprecision::mp_float& x);
+      bool near_one(const boost::multiprecision::mp_complex& z);
+
+      bool near_int(const double x);
+      bool near_int(const boost::multiprecision::mp_float& x);
+      bool near_int(const boost::multiprecision::mp_complex& z);
+
+      const boost::multiprecision::mp_float& two                     (void);
+      const boost::multiprecision::mp_float& three                   (void);
+      const boost::multiprecision::mp_float& four                    (void);
+      const boost::multiprecision::mp_float& five                    (void);
+      const boost::multiprecision::mp_float& six                     (void);
+      const boost::multiprecision::mp_float& seven                   (void);
+      const boost::multiprecision::mp_float& eight                   (void);
+      const boost::multiprecision::mp_float& nine                    (void);
+      const boost::multiprecision::mp_float& ten                     (void);
+      const boost::multiprecision::mp_float& twenty                  (void);
+      const boost::multiprecision::mp_float& thirty                  (void);
+      const boost::multiprecision::mp_float& forty                   (void);
+      const boost::multiprecision::mp_float& fifty                   (void);
+      const boost::multiprecision::mp_float& hundred                 (void);
+      const boost::multiprecision::mp_float& two_hundred             (void);
+      const boost::multiprecision::mp_float& three_hundred           (void);
+      const boost::multiprecision::mp_float& four_hundred            (void);
+      const boost::multiprecision::mp_float& five_hundred            (void);
+      const boost::multiprecision::mp_float& thousand                (void);
+      const boost::multiprecision::mp_float& two_k                   (void);
+      const boost::multiprecision::mp_float& three_k                 (void);
+      const boost::multiprecision::mp_float& four_k                  (void);
+      const boost::multiprecision::mp_float& five_k                  (void);
+      const boost::multiprecision::mp_float& ten_k                   (void);
+      const boost::multiprecision::mp_float& twenty_k                (void);
+      const boost::multiprecision::mp_float& thirty_k                (void);
+      const boost::multiprecision::mp_float& forty_k                 (void);
+      const boost::multiprecision::mp_float& fifty_k                 (void);
+      const boost::multiprecision::mp_float& hundred_k               (void);
+      const boost::multiprecision::mp_float& million                 (void);
+      const boost::multiprecision::mp_float& ten_M                   (void);
+      const boost::multiprecision::mp_float& hundred_M               (void);
+      const boost::multiprecision::mp_float& billion                 (void);
+      const boost::multiprecision::mp_float& trillion                (void);
+      const boost::multiprecision::mp_float& googol                  (void);
+      const boost::multiprecision::mp_float& int64_min               (void);
+      const boost::multiprecision::mp_float& int64_max               (void);
+      const boost::multiprecision::mp_float& int32_min               (void);
+      const boost::multiprecision::mp_float& int32_max               (void);
+      const boost::multiprecision::mp_float& unsigned_long_long_max  (void);
+      const boost::multiprecision::mp_float& signed_long_long_max    (void);
+      const boost::multiprecision::mp_float& signed_long_long_min    (void);
+      const boost::multiprecision::mp_float& double_min              (void);
+      const boost::multiprecision::mp_float& double_max              (void);
+      const boost::multiprecision::mp_float& long_double_min         (void);
+      const boost::multiprecision::mp_float& long_double_max         (void);
+      const boost::multiprecision::mp_float& one_minus               (void);
+      const boost::multiprecision::mp_float& tenth                   (void);
+      const boost::multiprecision::mp_float& eighth                  (void);
+      const boost::multiprecision::mp_float& sixteenth               (void);
+      const boost::multiprecision::mp_float& fifth                   (void);
+      const boost::multiprecision::mp_float& quarter                 (void);
+      const boost::multiprecision::mp_float& third                   (void);
+      const boost::multiprecision::mp_float& two_third               (void);
+      const boost::multiprecision::mp_float& four_third              (void);
+      const boost::multiprecision::mp_float& three_half              (void);
+      const boost::multiprecision::mp_float& sqrt2                   (void);
+      const boost::multiprecision::mp_float& sqrt3                   (void);
+      const boost::multiprecision::mp_float& pi                      (void);
+      const boost::multiprecision::mp_float& pi_half                 (void);
+      const boost::multiprecision::mp_float& pi_quarter              (void);
+      const boost::multiprecision::mp_float& pi_squared              (void);
+      const boost::multiprecision::mp_float& two_pi                  (void);
+      const boost::multiprecision::mp_float& sqrt_pi                 (void);
+      const boost::multiprecision::mp_float& degree                  (void);
+      const boost::multiprecision::mp_float& exp1                    (void);
+      const boost::multiprecision::mp_float& ln2                     (void);
+      const boost::multiprecision::mp_float& ln3                     (void);
+      const boost::multiprecision::mp_float& ln10                    (void);
+      const boost::multiprecision::mp_float& log10_2                 (void);
+      const boost::multiprecision::mp_float& golden_ratio            (void);
+      const boost::multiprecision::mp_float& euler_gamma             (void);
+      const boost::multiprecision::mp_float& catalan                 (void);
+      const boost::multiprecision::mp_float& khinchin                (void);
+      const boost::multiprecision::mp_float& glaisher                (void);
+      const boost::multiprecision::mp_float& extreme_value_skewness  (void);
+      const boost::multiprecision::mp_float& rayleigh_skewness       (void);
+      const boost::multiprecision::mp_float& rayleigh_kurtosis       (void);
+      const boost::multiprecision::mp_float& rayleigh_kurtosis_excess(void);
+
+      boost::multiprecision::mp_float pow2    (const boost::int64_t p);
+      boost::multiprecision::mp_float pown    (const boost::multiprecision::mp_float& x, const boost::int64_t p);
+      boost::multiprecision::mp_float inv     (const boost::multiprecision::mp_float& x);
+      boost::multiprecision::mp_float sqrt    (const boost::multiprecision::mp_float& x);
+      boost::multiprecision::mp_float cbrt    (const boost::multiprecision::mp_float& x);
+      boost::multiprecision::mp_float rootn   (const boost::multiprecision::mp_float& x, const boost::int32_t p);
+      boost::multiprecision::mp_float exp     (const boost::multiprecision::mp_float& x);
+      boost::multiprecision::mp_float log     (const boost::multiprecision::mp_float& x);
+      boost::multiprecision::mp_float log10   (const boost::multiprecision::mp_float& x);
+      boost::multiprecision::mp_float loga    (const boost::multiprecision::mp_float& a, const boost::multiprecision::mp_float& x);
+      boost::multiprecision::mp_float log1p   (const boost::multiprecision::mp_float& x);
+      boost::multiprecision::mp_float log1p1m2(const boost::multiprecision::mp_float& x);
+      boost::multiprecision::mp_float pow     (const boost::multiprecision::mp_float& x, const boost::multiprecision::mp_float& a);
+      void    sincos  (const boost::multiprecision::mp_float& x, boost::multiprecision::mp_float* const p_sin, boost::multiprecision::mp_float* const p_cos);
+      boost::multiprecision::mp_float sin     (const boost::multiprecision::mp_float& x);
+      boost::multiprecision::mp_float cos     (const boost::multiprecision::mp_float& x);
+      boost::multiprecision::mp_float tan     (const boost::multiprecision::mp_float& x);
+      boost::multiprecision::mp_float csc     (const boost::multiprecision::mp_float& x);
+      boost::multiprecision::mp_float sec     (const boost::multiprecision::mp_float& x);
+      boost::multiprecision::mp_float cot     (const boost::multiprecision::mp_float& x);
+      boost::multiprecision::mp_float asin    (const boost::multiprecision::mp_float& x);
+      boost::multiprecision::mp_float acos    (const boost::multiprecision::mp_float& x);
+      boost::multiprecision::mp_float atan    (const boost::multiprecision::mp_float& x);
+      boost::multiprecision::mp_float atan2   (const boost::multiprecision::mp_float& y, const boost::multiprecision::mp_float& x);
+      void    sinhcosh(const boost::multiprecision::mp_float& x, boost::multiprecision::mp_float* const p_sin, boost::multiprecision::mp_float* const p_cos);
+      boost::multiprecision::mp_float sinh    (const boost::multiprecision::mp_float& x);
+      boost::multiprecision::mp_float cosh    (const boost::multiprecision::mp_float& x);
+      boost::multiprecision::mp_float tanh    (const boost::multiprecision::mp_float& x);
+      boost::multiprecision::mp_float asinh   (const boost::multiprecision::mp_float& x);
+      boost::multiprecision::mp_float acosh   (const boost::multiprecision::mp_float& x);
+      boost::multiprecision::mp_float atanh   (const boost::multiprecision::mp_float& x);
+
+      boost::multiprecision::mp_float hyp0F0(const boost::multiprecision::mp_float& x);
+      boost::multiprecision::mp_float hyp0F1(const boost::multiprecision::mp_float& b, const boost::multiprecision::mp_float& x);
+      boost::multiprecision::mp_float hyp1F0(const boost::multiprecision::mp_float& a, const boost::multiprecision::mp_float& x);
+      boost::multiprecision::mp_float hyp1F1(const boost::multiprecision::mp_float& a, const boost::multiprecision::mp_float& b, const boost::multiprecision::mp_float& x);
+      boost::multiprecision::mp_float hyp2F0(const boost::multiprecision::mp_float& a, const boost::multiprecision::mp_float& b, const boost::multiprecision::mp_float& x);
+      boost::multiprecision::mp_float hyp2F1(const boost::multiprecision::mp_float& a, const boost::multiprecision::mp_float& b, const boost::multiprecision::mp_float& c, const boost::multiprecision::mp_float& x);
+      boost::multiprecision::mp_float hypPFQ(const std::deque<boost::multiprecision::mp_float>& a, const  std::deque<boost::multiprecision::mp_float>& b, const boost::multiprecision::mp_float& x);
+
+      boost::multiprecision::mp_float gamma          (const boost::multiprecision::mp_float& x);
+      boost::multiprecision::mp_float gamma_near_n   (const boost::int32_t n, const boost::multiprecision::mp_float& x);
+      boost::multiprecision::mp_float factorial      (const boost::uint32_t n);
+      boost::multiprecision::mp_float factorial2     (const  boost::int32_t n);
+      boost::multiprecision::mp_float binomial       (const boost::uint32_t n, const boost::uint32_t k);
+      boost::multiprecision::mp_float binomial       (const boost::uint32_t n, const boost::multiprecision::mp_float& y);
+      boost::multiprecision::mp_float binomial       (const boost::multiprecision::mp_float& x, const boost::uint32_t k);
+      boost::multiprecision::mp_float binomial       (const boost::multiprecision::mp_float& x, const boost::multiprecision::mp_float& y);
+      boost::multiprecision::mp_float pochhammer     (const boost::multiprecision::mp_float& x, const boost::uint32_t n);
+      boost::multiprecision::mp_float pochhammer     (const boost::multiprecision::mp_float& x, const boost::multiprecision::mp_float& a);
+      boost::multiprecision::mp_float bernoulli      (const boost::uint32_t n);
+      void    bernoulli_table(std::vector<boost::multiprecision::mp_float>& bn, const boost::uint32_t n);
+      void    prime          (const boost::uint32_t n, std::deque<boost::uint32_t>& primes);
+      boost::multiprecision::mp_float riemann_zeta   (const boost::int32_t n);
+      boost::multiprecision::mp_float riemann_zeta   (const boost::multiprecision::mp_float& s);
+    }
+  }
+
+  namespace boost
+  {
+    namespace multiprecision
+    {
+      inline boost::multiprecision::mp_float norm(const boost::multiprecision::mp_complex& z) { return z.norm(); }
+             boost::multiprecision::mp_float abs (const boost::multiprecision::mp_complex& z);
+             boost::multiprecision::mp_float arg (const boost::multiprecision::mp_complex& z);
+      inline boost::multiprecision::mp_float real(const boost::multiprecision::mp_complex& z) { return z.real(); }
+      inline boost::multiprecision::mp_float imag(const boost::multiprecision::mp_complex& z) { return z.imag(); }
+
+      inline boost::multiprecision::mp_complex conj(const boost::multiprecision::mp_complex& z) { return boost::multiprecision::mp_complex(z.real(), -z.imag()); }
+      inline boost::multiprecision::mp_complex iz  (const boost::multiprecision::mp_complex& z) { const boost::multiprecision::mp_float tmp(z.real()); return boost::multiprecision::mp_complex(-z.imag(), tmp); }
+
+      boost::multiprecision::mp_complex polar   (const boost::multiprecision::mp_float& mod, const boost::multiprecision::mp_float& arg);
+      boost::multiprecision::mp_complex sin     (const boost::multiprecision::mp_complex& z);
+      boost::multiprecision::mp_complex cos     (const boost::multiprecision::mp_complex& z);
+      boost::multiprecision::mp_complex tan     (const boost::multiprecision::mp_complex& z);
+      void       sincos  (const boost::multiprecision::mp_complex& z, boost::multiprecision::mp_complex* const p_sin, boost::multiprecision::mp_complex* const p_cos);
+      boost::multiprecision::mp_complex csc     (const boost::multiprecision::mp_complex& z);
+      boost::multiprecision::mp_complex sec     (const boost::multiprecision::mp_complex& z);
+      boost::multiprecision::mp_complex cot     (const boost::multiprecision::mp_complex& z);
+      boost::multiprecision::mp_complex asin    (const boost::multiprecision::mp_complex& z);
+      boost::multiprecision::mp_complex acos    (const boost::multiprecision::mp_complex& z);
+      boost::multiprecision::mp_complex atan    (const boost::multiprecision::mp_complex& z);
+      boost::multiprecision::mp_complex inv     (const boost::multiprecision::mp_complex& z);
+      boost::multiprecision::mp_complex sqrt    (const boost::multiprecision::mp_complex& z);
+      boost::multiprecision::mp_complex exp     (const boost::multiprecision::mp_complex& z);
+      boost::multiprecision::mp_complex log     (const boost::multiprecision::mp_complex& z);
+      boost::multiprecision::mp_complex log10   (const boost::multiprecision::mp_complex& z);
+      boost::multiprecision::mp_complex loga    (const boost::multiprecision::mp_complex& a, const boost::multiprecision::mp_complex& z);
+      boost::multiprecision::mp_complex pown    (const boost::multiprecision::mp_complex& z, const boost::int64_t p);
+      boost::multiprecision::mp_complex pow     (const boost::multiprecision::mp_complex& z, const boost::multiprecision::mp_complex& a);
+      boost::multiprecision::mp_complex rootn   (const boost::multiprecision::mp_complex& z, const boost::int32_t p);
+      boost::multiprecision::mp_complex sinh    (const boost::multiprecision::mp_complex& z);
+      boost::multiprecision::mp_complex cosh    (const boost::multiprecision::mp_complex& z);
+      boost::multiprecision::mp_complex tanh    (const boost::multiprecision::mp_complex& z);
+      void       sinhcosh(const boost::multiprecision::mp_complex& z, boost::multiprecision::mp_complex* const p_sinh, boost::multiprecision::mp_complex* const p_cosh);
+      boost::multiprecision::mp_complex asinh   (const boost::multiprecision::mp_complex& z);
+      boost::multiprecision::mp_complex acosh   (const boost::multiprecision::mp_complex& z);
+      boost::multiprecision::mp_complex atanh   (const boost::multiprecision::mp_complex& z);
+
+      boost::multiprecision::mp_complex gamma       (const boost::multiprecision::mp_complex& z);
+      boost::multiprecision::mp_complex pochhammer  (const boost::multiprecision::mp_complex& z, const boost::uint32_t n);
+      boost::multiprecision::mp_complex pochhammer  (const boost::multiprecision::mp_complex& z, const boost::multiprecision::mp_complex& a);
+      boost::multiprecision::mp_complex riemann_zeta(const boost::multiprecision::mp_complex& s);
+    }
+  }
+
+#endif // _MP_FLOAT_FUNCTIONS_HPP_

Added: sandbox/multiprecision/boost/mp_float_gmp.hpp
==============================================================================
--- (empty file)
+++ sandbox/multiprecision/boost/mp_float_gmp.hpp	2011-09-21 17:00:15 EDT (Wed, 21 Sep 2011)
@@ -0,0 +1,176 @@
+
+//          Copyright Christopher Kormanyos 2002 - 2011.
+// Distributed under the Boost Software License, Version 1.0.
+//    (See accompanying file LICENSE_1_0.txt or copy at
+//          http://www.boost.org/LICENSE_1_0.txt)
+
+// This work is based on an earlier work:
+// "Algorithm 910: A Portable C++ Multiple-Precision System for Special-Function Calculations",
+// in ACM TOMS, {VOL 37, ISSUE 4, (February 2011)} (C) ACM, 2011. http://doi.acm.org/10.1145/1916461.1916469
+
+#ifndef _MP_FLOAT_GMP_HPP_
+  #define _MP_FLOAT_GMP_HPP_
+
+  #if defined(__INTEL_COMPILER)
+    #pragma warning (disable:193)
+    #pragma warning (disable:981)
+  #endif
+
+  #if defined(_MSC_VER)
+    #pragma warning (disable:4127)
+  #endif
+
+  #include <cmath>
+  #include <climits>
+  #include <string>
+
+  #include <boost/multiprecision/mp_float_base.hpp>
+
+  // Declare the types of GMP.
+  extern "C"
+  {
+    typedef long int          mp_size_t;
+    typedef long int          mp_exp_t;
+    typedef unsigned long int mp_limb_t;
+
+    typedef struct struct__mpf_struct
+    {
+      int _mp_prec;
+      int _mp_size;
+      mp_exp_t _mp_exp;
+      mp_limb_t *_mp_d;
+    }
+    __mpf_struct;
+
+    typedef       __mpf_struct mpf_t[1];
+    typedef       __mpf_struct* mpf_ptr;
+    typedef const __mpf_struct* mpf_srcptr;
+  }
+
+  namespace boost
+  {
+    namespace multiprecision
+    {
+      class mp_float_gmp : public mp_float_base
+      {
+      public:
+        static const boost::int32_t mp_digits    = static_cast<boost::int32_t>((static_cast<signed long long>(mp_digits10) * 2136LL) / 643LL);
+        static const boost::int32_t mp_radix     = static_cast<boost::int32_t>(2);
+
+        static const boost::int64_t mp_max_exp   = static_cast<boost::int64_t>(LONG_MAX - 31LL); // TBD: Ensure that (boost::int64_t >= long)
+        static const boost::int64_t mp_min_exp   = static_cast<boost::int64_t>(LONG_MIN + 31LL); // TBD: Ensure that (boost::int64_t >= long)
+        static const boost::int64_t mp_max_exp10 = static_cast<boost::int64_t>((static_cast<signed long long>(mp_max_exp) * 643LL) / 2136LL);
+        static const boost::int64_t mp_min_exp10 = static_cast<boost::int64_t>((static_cast<signed long long>(mp_min_exp) * 643LL) / 2136LL);
+
+      private:
+        static const boost::int32_t mp_digits2 = static_cast<boost::int32_t>((static_cast<signed long long>(mp_max_digits10) * 2136LL) / 643LL);
+
+        typedef enum enum_fpclass
+        {
+          mp_finite,
+          mp_inf_pos,
+          mp_inf_neg,
+          mp_NaN
+        }
+        t_fpclass;
+
+        t_fpclass fpclass;
+        boost::int32_t     prec_elem;
+        ::mpf_t   rop;
+
+      public:
+        // Constructors
+        mp_float_gmp();
+        mp_float_gmp(const char n);
+        mp_float_gmp(const signed char n);
+        mp_float_gmp(const unsigned char n);
+        mp_float_gmp(const wchar_t n);
+        mp_float_gmp(const signed short n);
+        mp_float_gmp(const unsigned short n);
+        mp_float_gmp(const signed int n);
+        mp_float_gmp(const unsigned int n);
+        mp_float_gmp(const signed long n);
+        mp_float_gmp(const unsigned long n);
+        mp_float_gmp(const signed long long n);
+        mp_float_gmp(const unsigned long long n);
+        mp_float_gmp(const float f);
+        mp_float_gmp(const double d);
+        mp_float_gmp(const long double ld);
+        mp_float_gmp(const char* const s);
+        mp_float_gmp(const std::string& str);
+
+        mp_float_gmp(const mp_float_gmp& f);
+        mp_float_gmp(const double mantissa, const boost::int64_t exponent);
+
+        virtual ~mp_float_gmp();
+
+      private:
+        explicit mp_float_gmp(const ::mpf_t& op);
+
+      public:
+        virtual boost::int32_t cmp(const mp_float_gmp& v) const;
+
+        // Specific special values.
+        virtual const mp_float_gmp& my_value_nan(void) const;
+        virtual const mp_float_gmp& my_value_inf(void) const;
+
+        virtual void precision(const boost::int32_t prec_digits);
+
+        // Basic operations.
+        virtual mp_float_gmp& operator= (const mp_float_gmp& v);
+        virtual mp_float_gmp& operator+=(const mp_float_gmp& v);
+        virtual mp_float_gmp& operator-=(const mp_float_gmp& v);
+        virtual mp_float_gmp& operator*=(const mp_float_gmp& v);
+        virtual mp_float_gmp& operator/=(const mp_float_gmp& v);
+        virtual mp_float_gmp& add_unsigned_long_long(const unsigned long long n);
+        virtual mp_float_gmp& sub_unsigned_long_long(const unsigned long long n);
+        virtual mp_float_gmp& mul_unsigned_long_long(const unsigned long long n);
+        virtual mp_float_gmp& div_unsigned_long_long(const unsigned long long n);
+
+        // Elementary primitives.
+        virtual mp_float_gmp& calculate_inv (void);
+        virtual mp_float_gmp& calculate_sqrt(void);
+        virtual mp_float_gmp& negate(void);
+
+        // Comparison functions
+        virtual bool isnan   (void) const { return  (fpclass == mp_NaN); }
+        virtual bool isinf   (void) const { return ((fpclass == mp_inf_pos) || (fpclass == mp_inf_neg)); }
+        virtual bool isfinite(void) const { return  (fpclass == mp_finite); }
+
+        virtual bool iszero (void) const;
+        virtual bool isone  (void) const;
+        virtual bool isint  (void) const;
+        virtual bool isneg  (void) const;
+
+        // Operators pre-increment and pre-decrement
+        virtual mp_float_gmp& operator++(void);
+        virtual mp_float_gmp& operator--(void);
+
+        // Conversion routines
+        virtual void               extract_parts             (double& mantissa, boost::int64_t& exponent) const;
+        virtual double             extract_double            (void) const;
+        virtual long double        extract_long_double       (void) const;
+        virtual signed long long   extract_signed_long_long  (void) const;
+        virtual unsigned long long extract_unsigned_long_long(void) const;
+        virtual mp_float_gmp            extract_integer_part      (void) const;
+        virtual mp_float_gmp            extract_decimal_part      (void) const;
+
+      private:
+        static const boost::int64_t& max_exp2(void);
+        static const boost::int64_t& min_exp2(void);
+
+        static void init(void);
+
+        void from_unsigned_long_long(const unsigned long long u);
+        void from_unsigned_long(const unsigned long u);
+
+        virtual bool rd_string(const char* const s);
+
+        virtual boost::int64_t get_order_exact(void) const;
+        virtual boost::int64_t get_order_fast(void) const;
+        virtual void get_output_string(std::string& str, boost::int64_t& my_exp, const std::size_t number_of_digits) const;
+      };
+    }
+  }
+
+#endif // _MP_FLOAT_GMP_HPP_

Added: sandbox/multiprecision/boost/mp_float_mpfr.hpp
==============================================================================
--- (empty file)
+++ sandbox/multiprecision/boost/mp_float_mpfr.hpp	2011-09-21 17:00:15 EDT (Wed, 21 Sep 2011)
@@ -0,0 +1,207 @@
+
+//          Copyright Christopher Kormanyos 2002 - 2011.
+// Distributed under the Boost Software License, Version 1.0.
+//    (See accompanying file LICENSE_1_0.txt or copy at
+//          http://www.boost.org/LICENSE_1_0.txt)
+
+// This work is based on an earlier work:
+// "Algorithm 910: A Portable C++ Multiple-Precision System for Special-Function Calculations",
+// in ACM TOMS, {VOL 37, ISSUE 4, (February 2011)} (C) ACM, 2011. http://doi.acm.org/10.1145/1916461.1916469
+
+#ifndef _MP_FLOAT_MPFR_HPP_
+  #define _MP_FLOAT_MPFR_HPP_
+
+  #include <cmath>
+  #include <string>
+  #include <climits>
+
+  #include <boost/multiprecision/mp_float_base.hpp>
+
+  // Declare the types of MPFR.
+  extern "C"
+  {
+    typedef unsigned long     mpfr_prec_t;
+    typedef int               mpfr_sign_t;
+    typedef long int          mp_exp_t;
+    typedef unsigned long int mp_limb_t;
+
+    #define mp_prec_t mpfr_prec_t
+    #define mp_rnd_t  mpfr_rnd_t
+
+    typedef struct
+    {
+      mpfr_prec_t _mpfr_prec;
+      mpfr_sign_t _mpfr_sign;
+      mp_exp_t    _mpfr_exp;
+      mp_limb_t*  _mpfr_d;
+    }
+    __mpfr_struct;
+
+    #define __gmp_const const
+
+    typedef __mpfr_struct mpfr_t[1];
+    typedef __mpfr_struct* mpfr_ptr;
+    typedef __gmp_const __mpfr_struct* mpfr_srcptr;
+
+    typedef enum
+    {
+      GMP_RNDNA   = -1,
+      GMP_RNDN    =  0,
+      GMP_RNDZ    =  1,
+      GMP_RNDU    =  2,
+      GMP_RNDD    =  3,
+      GMP_RND_MAX =  4
+    }
+    mpfr_rnd_t;
+  }
+
+  namespace boost
+  {
+    namespace multiprecision
+    {
+      class mp_float_mpfr : public mp_float_base
+      {
+      public:
+        static const boost::int32_t mp_digits    = static_cast<boost::int32_t>((static_cast<signed long long>(mp_digits10) * 2136LL) / 643LL);
+        static const boost::int32_t mp_radix     = 2;
+
+        static const boost::int64_t mp_max_exp   = static_cast<boost::int64_t>(LONG_MAX / static_cast<signed long>(2L)); // TBD: Ensure that (boost::int64_t >= long)
+        static const boost::int64_t mp_min_exp   = static_cast<boost::int64_t>(LONG_MIN / static_cast<signed long>(2L)); // TBD: Ensure that (boost::int64_t >= long)
+        static const boost::int64_t mp_max_exp10 = static_cast<boost::int64_t>((static_cast<signed long long>(mp_max_exp) * 643LL) / 2136LL);
+        static const boost::int64_t mp_min_exp10 = static_cast<boost::int64_t>((static_cast<signed long long>(mp_min_exp) * 643LL) / 2136LL);
+
+      private:
+        static const boost::int32_t mp_digits2 = static_cast<boost::int32_t>((static_cast<signed long long>(mp_max_digits10) * 2136LL) / 643LL);
+        ::mpfr_t rop;
+
+      public:
+        // Constructors
+        mp_float_mpfr();
+        mp_float_mpfr(const char n);
+        mp_float_mpfr(const signed char n);
+        mp_float_mpfr(const unsigned char n);
+        mp_float_mpfr(const wchar_t n);
+        mp_float_mpfr(const signed short n);
+        mp_float_mpfr(const unsigned short n);
+        mp_float_mpfr(const signed int n);
+        mp_float_mpfr(const unsigned int n);
+        mp_float_mpfr(const signed long n);
+        mp_float_mpfr(const unsigned long n);
+        mp_float_mpfr(const signed long long n);
+        mp_float_mpfr(const unsigned long long n);
+        mp_float_mpfr(const float f);
+        mp_float_mpfr(const double d);
+        mp_float_mpfr(const long double ld);
+        mp_float_mpfr(const char* const s);
+        mp_float_mpfr(const std::string& str);
+
+        mp_float_mpfr(const mp_float_mpfr& f);
+        mp_float_mpfr(const double mantissa, const boost::int64_t exponent);
+
+        virtual ~mp_float_mpfr();
+
+        virtual boost::int32_t cmp(const mp_float_mpfr& v) const;
+
+        // Specific special values.
+        virtual const mp_float_mpfr& my_value_nan(void) const;
+        virtual const mp_float_mpfr& my_value_inf(void) const;
+
+        virtual void precision(const boost::int32_t) { }
+
+        // Basic operations.
+        virtual mp_float_mpfr& operator= (const mp_float_mpfr& v);
+        virtual mp_float_mpfr& operator+=(const mp_float_mpfr& v);
+        virtual mp_float_mpfr& operator-=(const mp_float_mpfr& v);
+        virtual mp_float_mpfr& operator*=(const mp_float_mpfr& v);
+        virtual mp_float_mpfr& operator/=(const mp_float_mpfr& v);
+        virtual mp_float_mpfr& add_unsigned_long_long(const unsigned long long n);
+        virtual mp_float_mpfr& sub_unsigned_long_long(const unsigned long long n);
+        virtual mp_float_mpfr& mul_unsigned_long_long(const unsigned long long n);
+        virtual mp_float_mpfr& div_unsigned_long_long(const unsigned long long n);
+
+        // Elementary primitives.
+        virtual mp_float_mpfr& calculate_inv (void);
+        virtual mp_float_mpfr& calculate_sqrt(void);
+        virtual mp_float_mpfr& negate(void);
+
+        // Comparison functions
+        virtual bool isnan   (void) const;
+        virtual bool isinf   (void) const;
+        virtual bool isfinite(void) const;
+
+        virtual bool iszero  (void) const;
+        virtual bool isone   (void) const;
+        virtual bool isint   (void) const;
+        virtual bool isneg   (void) const;
+
+        // Operators pre-increment and pre-decrement.
+        virtual mp_float_mpfr& operator++(void);
+        virtual mp_float_mpfr& operator--(void);
+
+        // Conversion routines.
+        virtual void               extract_parts             (double& mantissa, boost::int64_t& exponent) const;
+        virtual double             extract_double            (void) const;
+        virtual long double        extract_long_double       (void) const;
+        virtual signed long long   extract_signed_long_long  (void) const;
+        virtual unsigned long long extract_unsigned_long_long(void) const;
+        virtual mp_float_mpfr            extract_integer_part      (void) const;
+        virtual mp_float_mpfr            extract_decimal_part      (void) const;
+
+        static mp_float_mpfr my_cbrt         (const mp_float_mpfr& x);
+        static mp_float_mpfr my_rootn        (const mp_float_mpfr& x, const boost::uint32_t p);
+        static mp_float_mpfr my_exp          (const mp_float_mpfr& x);
+        static mp_float_mpfr my_log          (const mp_float_mpfr& x);
+        static mp_float_mpfr my_sin          (const mp_float_mpfr& x);
+        static mp_float_mpfr my_cos          (const mp_float_mpfr& x);
+        static mp_float_mpfr my_tan          (const mp_float_mpfr& x);
+        static mp_float_mpfr my_asin         (const mp_float_mpfr& x);
+        static mp_float_mpfr my_acos         (const mp_float_mpfr& x);
+        static mp_float_mpfr my_atan         (const mp_float_mpfr& x);
+        static mp_float_mpfr my_sinh         (const mp_float_mpfr& x);
+        static mp_float_mpfr my_cosh         (const mp_float_mpfr& x);
+        static mp_float_mpfr my_tanh         (const mp_float_mpfr& x);
+        static mp_float_mpfr my_asinh        (const mp_float_mpfr& x);
+        static mp_float_mpfr my_acosh        (const mp_float_mpfr& x);
+        static mp_float_mpfr my_atanh        (const mp_float_mpfr& x);
+        static mp_float_mpfr my_gamma        (const mp_float_mpfr& x);
+        static mp_float_mpfr my_riemann_zeta (const mp_float_mpfr& x);
+        static mp_float_mpfr my_cyl_bessel_jn(const boost::int32_t n, const mp_float_mpfr& x);
+        static mp_float_mpfr my_cyl_bessel_yn(const boost::int32_t n, const mp_float_mpfr& x);
+
+        virtual bool has_its_own_cbrt         (void) const { return true; }
+        virtual bool has_its_own_rootn        (void) const { return true; }
+        virtual bool has_its_own_exp          (void) const { return true; }
+        virtual bool has_its_own_log          (void) const { return true; }
+        virtual bool has_its_own_sin          (void) const { return true; }
+        virtual bool has_its_own_cos          (void) const { return true; }
+        virtual bool has_its_own_tan          (void) const { return true; }
+        virtual bool has_its_own_asin         (void) const { return true; }
+        virtual bool has_its_own_acos         (void) const { return true; }
+        virtual bool has_its_own_atan         (void) const { return true; }
+        virtual bool has_its_own_sinh         (void) const { return true; }
+        virtual bool has_its_own_cosh         (void) const { return true; }
+        virtual bool has_its_own_tanh         (void) const { return true; }
+        virtual bool has_its_own_asinh        (void) const { return true; }
+        virtual bool has_its_own_acosh        (void) const { return true; }
+        virtual bool has_its_own_atanh        (void) const { return true; }
+        virtual bool has_its_own_gamma        (void) const { return false; }
+        virtual bool has_its_own_riemann_zeta (void) const { return false; }
+        virtual bool has_its_own_cyl_bessel_jn(void) const { return false; }
+        virtual bool has_its_own_cyl_bessel_yn(void) const { return false; }
+
+      private:
+        static void init(void);
+
+        void from_unsigned_long_long(const unsigned long long u);
+        void from_unsigned_long(const unsigned long u);
+
+        virtual bool rd_string(const char* const s);
+
+        virtual boost::int64_t get_order_exact(void) const;
+        virtual boost::int64_t get_order_fast(void) const;
+        virtual void get_output_string(std::string& str, boost::int64_t& my_exp, const std::size_t number_of_digits) const;
+      };
+    }
+  }
+
+#endif // _MP_FLOAT_MPFR_HPP_

Added: sandbox/multiprecision/boost/utility/util_alternating_sum.hpp
==============================================================================
--- (empty file)
+++ sandbox/multiprecision/boost/utility/util_alternating_sum.hpp	2011-09-21 17:00:15 EDT (Wed, 21 Sep 2011)
@@ -0,0 +1,45 @@
+
+//          Copyright Christopher Kormanyos 2002 - 2011.
+// Distributed under the Boost Software License, Version 1.0.
+//    (See accompanying file LICENSE_1_0.txt or copy at
+//          http://www.boost.org/LICENSE_1_0.txt)
+
+// This work is based on an earlier work:
+// "Algorithm 910: A Portable C++ Multiple-Precision System for Special-Function Calculations",
+// in ACM TOMS, {VOL 37, ISSUE 4, (February 2011)} (C) ACM, 2011. http://doi.acm.org/10.1145/1916461.1916469
+
+#ifndef _UTIL_ALTERNATING_SUM_2010_01_11_HPP_
+  #define _UTIL_ALTERNATING_SUM_2010_01_11_HPP_
+
+  namespace boost
+  {
+    namespace multiprecision
+    {
+      namespace utility
+      {
+        template<typename T1, typename T2 = T1> struct alternating_sum
+        {
+        private:
+
+          const alternating_sum& operator=(const alternating_sum&);
+
+          bool b_neg_term;
+          const T2 initial;
+
+        public:
+
+          alternating_sum(const bool b_neg = false, const T2& init = T2(0)) : b_neg_term(b_neg),
+                                                                              initial   (init) { }
+                                                                          
+          T1 operator()(const T1& sum, const T2& ck)
+          {
+            const T1 the_sum = (!b_neg_term ? (sum + ck) : (sum - ck));
+            b_neg_term = !b_neg_term;
+            return the_sum + initial;
+          }
+        };
+      }
+    }
+  }
+
+#endif // _UTIL_ALTERNATING_SUM_2010_01_11_HPP_

Added: sandbox/multiprecision/boost/utility/util_coefficient_expansion.hpp
==============================================================================
--- (empty file)
+++ sandbox/multiprecision/boost/utility/util_coefficient_expansion.hpp	2011-09-21 17:00:15 EDT (Wed, 21 Sep 2011)
@@ -0,0 +1,55 @@
+
+//          Copyright Christopher Kormanyos 2002 - 2011.
+// Distributed under the Boost Software License, Version 1.0.
+//    (See accompanying file LICENSE_1_0.txt or copy at
+//          http://www.boost.org/LICENSE_1_0.txt)
+
+// This work is based on an earlier work:
+// "Algorithm 910: A Portable C++ Multiple-Precision System for Special-Function Calculations",
+// in ACM TOMS, {VOL 37, ISSUE 4, (February 2011)} (C) ACM, 2011. http://doi.acm.org/10.1145/1916461.1916469
+
+#ifndef _UTIL_COEFFICIENT_EXPANSION_2009_11_23_HPP_
+  #define _UTIL_COEFFICIENT_EXPANSION_2009_11_23_HPP_
+
+  namespace boost
+  {
+    namespace multiprecision
+    {
+      namespace utility
+      {
+        template<typename T1, typename T2 = T1> struct coefficient_expansion
+        {
+        private:
+
+          const T1 x_expand;
+          T1       x_expand_pow_k;
+
+        private:
+
+          const coefficient_expansion& operator=(const coefficient_expansion&);
+
+          static const T2& one_t2(void)
+          {
+            static const T2 val_t2(1);
+            return val_t2;
+          }
+
+        public:
+
+          coefficient_expansion(const T1& expand, const T2& init = one_t2()) : x_expand      (expand),
+                                                                               x_expand_pow_k(init) { }
+
+          T1 operator()(const T1& sum, const T2& ck)
+          {
+            const T1 ck_x_pow_k(ck * x_expand_pow_k);
+
+            x_expand_pow_k *= x_expand;
+
+            return sum + ck_x_pow_k;
+          }
+        };
+      }
+    }
+  }
+
+#endif // _UTIL_COEFFICIENT_EXPANSION_2009_11_23_HPP_

Added: sandbox/multiprecision/boost/utility/util_digit_scale.hpp
==============================================================================
--- (empty file)
+++ sandbox/multiprecision/boost/utility/util_digit_scale.hpp	2011-09-21 17:00:15 EDT (Wed, 21 Sep 2011)
@@ -0,0 +1,25 @@
+
+//          Copyright Christopher Kormanyos 2002 - 2011.
+// Distributed under the Boost Software License, Version 1.0.
+//    (See accompanying file LICENSE_1_0.txt or copy at
+//          http://www.boost.org/LICENSE_1_0.txt)
+
+// This work is based on an earlier work:
+// "Algorithm 910: A Portable C++ Multiple-Precision System for Special-Function Calculations",
+// in ACM TOMS, {VOL 37, ISSUE 4, (February 2011)} (C) ACM, 2011. http://doi.acm.org/10.1145/1916461.1916469
+
+#ifndef _UTIL_DIGIT_SCALE_2009_11_01_HPP_
+  #define _UTIL_DIGIT_SCALE_2009_11_01_HPP_
+
+  namespace boost
+  {
+    namespace multiprecision
+    {
+      namespace utility
+      {
+        const double& digit_scale(void);
+      }
+    }
+  }
+
+#endif // _UTIL_POW_2009_01_25_HPP_

Added: sandbox/multiprecision/boost/utility/util_find_root_base.hpp
==============================================================================
--- (empty file)
+++ sandbox/multiprecision/boost/utility/util_find_root_base.hpp	2011-09-21 17:00:15 EDT (Wed, 21 Sep 2011)
@@ -0,0 +1,38 @@
+
+//          Copyright Christopher Kormanyos 2002 - 2011.
+// Distributed under the Boost Software License, Version 1.0.
+//    (See accompanying file LICENSE_1_0.txt or copy at
+//          http://www.boost.org/LICENSE_1_0.txt)
+
+// This work is based on an earlier work:
+// "Algorithm 910: A Portable C++ Multiple-Precision System for Special-Function Calculations",
+// in ACM TOMS, {VOL 37, ISSUE 4, (February 2011)} (C) ACM, 2011. http://doi.acm.org/10.1145/1916461.1916469
+
+#ifndef _UTIL_FIND_ROOT_BASE_2009_10_31_HPP_
+  #define _UTIL_FIND_ROOT_BASE_2009_10_31_HPP_
+
+  #include "util_ranged_function_operation.hpp"
+
+  namespace boost
+  {
+    namespace multiprecision
+    {
+      namespace utility
+      {
+        template<typename T> class find_root_base : public ranged_function_operation<T>
+        {
+        protected:
+    
+          find_root_base(const T& lo,
+                       const T& hi,
+                       const T& tol) : ranged_function_operation<T>(lo, hi, tol) { }
+
+        public:
+
+          virtual ~find_root_base() { }
+        };
+      }
+    }
+  }
+
+#endif // _UTIL_FIND_ROOT_BASE_2009_10_31_HPP_

Added: sandbox/multiprecision/boost/utility/util_find_root_bisect.hpp
==============================================================================
--- (empty file)
+++ sandbox/multiprecision/boost/utility/util_find_root_bisect.hpp	2011-09-21 17:00:15 EDT (Wed, 21 Sep 2011)
@@ -0,0 +1,106 @@
+
+//          Copyright Christopher Kormanyos 2002 - 2011.
+// Distributed under the Boost Software License, Version 1.0.
+//    (See accompanying file LICENSE_1_0.txt or copy at
+//          http://www.boost.org/LICENSE_1_0.txt)
+
+// This work is based on an earlier work:
+// "Algorithm 910: A Portable C++ Multiple-Precision System for Special-Function Calculations",
+// in ACM TOMS, {VOL 37, ISSUE 4, (February 2011)} (C) ACM, 2011. http://doi.acm.org/10.1145/1916461.1916469
+
+#ifndef _UTIL_FIND_ROOT_BISECT_2009_10_31_HPP_
+  #define _UTIL_FIND_ROOT_BISECT_2009_10_31_HPP_
+
+  #include <utility/util_find_root_base.h>
+  #include <boost/multiprecision/mp_float.hpp>
+
+  namespace boost
+  {
+    namespace multiprecision
+    {
+      namespace utility
+      {
+        template<typename T> class find_root_bisect : public find_root_base<T>
+        {
+        protected:
+    
+          find_root_bisect(const T& lo,
+                         const T& hi,
+                         const T& tol) : find_root_base<T>(lo, hi, tol) { }
+
+        public:
+
+          virtual ~find_root_bisect() { }
+
+        private:
+
+          virtual T my_operation(void) const
+          {
+            // Bisection method as described in Numerical Recipes in C++ 2nd Ed., chapter 9.1.
+            // The program on page 358 was taken directly from the book and slightly modified
+            // to improve adherence with standard C++ coding practices.
+
+            function_operation<T>::op_ok = true;
+
+            T lo = ranged_function_operation<T>::xlo;
+            T hi = ranged_function_operation<T>::xhi;
+
+            const T f = function(lo);
+
+            static const T t_zero = T(0);
+
+            // Make sure that there is at least one root in the interval.
+            if(f * function(ranged_function_operation<T>::xhi) >= t_zero)
+            {
+              return t_zero;
+            }
+
+            // Orient the search such that f > 0 lies at x + dx.
+            T dx;
+            T rt;
+
+            if(f < t_zero)
+            {
+              dx = hi - lo;
+              rt = lo;
+            }
+            else
+            {
+              dx = lo - hi;
+              rt = hi;
+            }
+
+            // Bisection iteration loop, maximum 2048 times.
+            for(boost::uint32_t i = static_cast<boost::uint32_t>(0u); i < static_cast<boost::uint32_t>(2048u); i++)
+            {
+              dx /= static_cast<boost::int32_t>(2);
+
+              const T x_mid = rt + dx;
+              const T f_mid = function(x_mid);
+
+              if(f_mid <= t_zero)
+              {
+                rt = x_mid;
+              }
+
+              // Check for convergence to within a tolerance.
+              const T dx_abs = ((dx < t_zero) ? -dx : dx);
+
+              if(dx_abs < ranged_function_operation<T>::eps || boost::multiprecision::iszero(f_mid))
+              {
+                // Return root.
+                return rt;
+              }
+            }
+
+            // Bisection iteration did not converge.
+            function_operation<T>::op_ok = false;
+
+            return t_zero;
+          }
+        };
+      }
+    }
+  }
+
+#endif // _UTIL_FIND_ROOT_BISECT_2009_10_31_HPP_

Added: sandbox/multiprecision/boost/utility/util_find_root_newton_raphson.hpp
==============================================================================
--- (empty file)
+++ sandbox/multiprecision/boost/utility/util_find_root_newton_raphson.hpp	2011-09-21 17:00:15 EDT (Wed, 21 Sep 2011)
@@ -0,0 +1,78 @@
+
+//          Copyright Christopher Kormanyos 2002 - 2011.
+// Distributed under the Boost Software License, Version 1.0.
+//    (See accompanying file LICENSE_1_0.txt or copy at
+//          http://www.boost.org/LICENSE_1_0.txt)
+
+// This work is based on an earlier work:
+// "Algorithm 910: A Portable C++ Multiple-Precision System for Special-Function Calculations",
+// in ACM TOMS, {VOL 37, ISSUE 4, (February 2011)} (C) ACM, 2011. http://doi.acm.org/10.1145/1916461.1916469
+
+#ifndef _UTIL_FIND_ROOT_NEWTON_RAPHSON_2009_10_27_HPP_
+  #define _UTIL_FIND_ROOT_NEWTON_RAPHSON_2009_10_27_HPP_
+
+  #include "util_find_root_base.hpp"
+
+  namespace boost
+  {
+    namespace multiprecision
+    {
+      namespace utility
+      {
+        template<typename T> class find_root_newton_raphson : public find_root_base<T>
+        {
+        protected:
+    
+          find_root_newton_raphson(const T& lo,
+                                const T& hi,
+                                const T& tol) : find_root_base<T>(lo, hi, tol) { }
+
+        public:
+
+          virtual ~find_root_newton_raphson() { }
+
+        public:
+
+          void function_derivative(const T& x, T& f, T& d) const { my_function_derivative(x, f, d); }
+
+        private:
+
+          virtual void my_function_derivative(const T& x, T& f, T& d) const = 0;
+
+          virtual T my_operation(void) const
+          {
+            function_operation<T>::op_ok = true;
+
+            T x = (find_root_base<T>::xlo + find_root_base<T>::xhi) / static_cast<boost::int32_t>(2);
+
+            for(boost::int32_t j = static_cast<boost::int32_t>(0); j < static_cast<boost::int32_t>(256); j++)
+            {
+              T f, d;
+              my_function_derivative(x, f, d);
+
+              const T dx = f / d;
+
+              x -= dx;
+
+              if((find_root_base<T>::xlo - x) * (x - find_root_base<T>::xhi) < T(0))
+              {
+                function_operation<T>::op_ok = false;
+                break;
+              }
+
+              const T delta = ((dx < T(0)) ? -dx : dx);
+
+              if(delta < find_root_base<T>::eps)
+              {
+                break;
+              }
+            }
+
+            return (function_operation<T>::op_ok ? x : T(0));
+          }
+        };
+      }
+    }
+  }
+
+#endif // _UTIL_FIND_ROOT_NEWTON_RAPHSON_2009_10_27_HPP_

Added: sandbox/multiprecision/boost/utility/util_function_base.hpp
==============================================================================
--- (empty file)
+++ sandbox/multiprecision/boost/utility/util_function_base.hpp	2011-09-21 17:00:15 EDT (Wed, 21 Sep 2011)
@@ -0,0 +1,40 @@
+
+//          Copyright Christopher Kormanyos 2002 - 2011.
+// Distributed under the Boost Software License, Version 1.0.
+//    (See accompanying file LICENSE_1_0.txt or copy at
+//          http://www.boost.org/LICENSE_1_0.txt)
+
+// This work is based on an earlier work:
+// "Algorithm 910: A Portable C++ Multiple-Precision System for Special-Function Calculations",
+// in ACM TOMS, {VOL 37, ISSUE 4, (February 2011)} (C) ACM, 2011. http://doi.acm.org/10.1145/1916461.1916469
+
+#ifndef _UTIL_FUNCTION_BASE_2009_10_27_HPP_
+  #define _UTIL_FUNCTION_BASE_2009_10_27_HPP_
+
+  namespace boost
+  {
+    namespace multiprecision
+    {
+      namespace utility
+      {
+        template<typename T> class function_base
+        {
+        public:
+          virtual ~function_base() { }
+
+          T function(const T& x) const { return my_function(x); }
+
+        protected:
+          function_base() { }
+
+        private:
+          const function_base& operator=(const function_base&);
+          function_base(const function_base&);
+
+          virtual T my_function(const T&) const = 0;
+        };
+      }
+    }
+  }
+
+#endif // _UTIL_FUNCTION_BASE_2009_10_27_HPP_

Added: sandbox/multiprecision/boost/utility/util_function_derivative.hpp
==============================================================================
--- (empty file)
+++ sandbox/multiprecision/boost/utility/util_function_derivative.hpp	2011-09-21 17:00:15 EDT (Wed, 21 Sep 2011)
@@ -0,0 +1,96 @@
+
+//          Copyright Christopher Kormanyos 2002 - 2011.
+// Distributed under the Boost Software License, Version 1.0.
+//    (See accompanying file LICENSE_1_0.txt or copy at
+//          http://www.boost.org/LICENSE_1_0.txt)
+
+// This work is based on an earlier work:
+// "Algorithm 910: A Portable C++ Multiple-Precision System for Special-Function Calculations",
+// in ACM TOMS, {VOL 37, ISSUE 4, (February 2011)} (C) ACM, 2011. http://doi.acm.org/10.1145/1916461.1916469
+
+#ifndef _UTIL_FUNCTION_DERIVATIVE_2009_11_17_HPP_
+  #define _UTIL_FUNCTION_DERIVATIVE_2009_11_17_HPP_
+
+  #include <string>
+  #include <sstream>
+
+  #include "util_function_operation.h"
+  #include <boost/lexical_cast.hpp>
+
+  namespace boost
+  {
+    namespace multiprecision
+    {
+      namespace utility
+      {
+        template<typename T> class function_derivative : public function_operation<T>
+        {
+        private:
+
+          static const T& my_tol(void)
+          {
+            static bool is_init = false;
+
+            static T val_tol;
+
+            if(!is_init)
+            {
+              is_init = true;
+
+              // Set the default tolerance to be approximately 10^[-(digits10 * 1.15)/5].
+              static const double      tx = (static_cast<double>(std::numeric_limits<T>::digits10) * 1.15) / 5.0;
+              static const std::size_t tn = static_cast<std::size_t>(tx + 0.5);
+
+              std::stringstream ss;
+
+              ss << "1E-" + boost::lexical_cast(tn);
+
+              ss >> val_tol;
+            }
+
+            static const T the_tol = val_tol;
+
+            return the_tol;
+          }
+
+        protected:
+
+          const T my_x;
+          const T my_dx;
+
+        protected:
+
+          function_derivative(const T& x, const T& dx = my_tol()) : my_x(x), my_dx(dx) { }
+
+        public:
+
+          virtual ~function_derivative() { }
+
+        protected:
+
+          virtual T my_operation(void) const
+          {
+            function_operation<T>::op_ok = true;
+
+            // Compute the function derivative using a three point
+            // central difference rule of O(dx^6).
+            const T dx1 = my_dx;
+            const T dx2 = dx1 * static_cast<boost::int32_t>(2);
+            const T dx3 = dx1 * static_cast<boost::int32_t>(3);
+
+            const T m1 = (function(my_x + dx1) - function(my_x - dx1)) / static_cast<boost::int32_t>(2);
+            const T m2 = (function(my_x + dx2) - function(my_x - dx2)) / static_cast<boost::int32_t>(4);
+            const T m3 = (function(my_x + dx3) - function(my_x - dx3)) / static_cast<boost::int32_t>(6);
+
+            const T fifteen_m1 = static_cast<boost::int32_t>(15) * m1;
+            const T six_m2     = static_cast<boost::int32_t>( 6) * m2;
+            const T ten_dx1    = static_cast<boost::int32_t>(10) * dx1;
+
+            return ((fifteen_m1 - six_m2) + m3) / ten_dx1;
+          }
+        };
+      }
+    }
+  }
+
+#endif // _UTIL_FUNCTION_DERIVATIVE_2009_11_17_HPP_

Added: sandbox/multiprecision/boost/utility/util_function_operation.hpp
==============================================================================
--- (empty file)
+++ sandbox/multiprecision/boost/utility/util_function_operation.hpp	2011-09-21 17:00:15 EDT (Wed, 21 Sep 2011)
@@ -0,0 +1,47 @@
+
+//          Copyright Christopher Kormanyos 2002 - 2011.
+// Distributed under the Boost Software License, Version 1.0.
+//    (See accompanying file LICENSE_1_0.txt or copy at
+//          http://www.boost.org/LICENSE_1_0.txt)
+
+// This work is based on an earlier work:
+// "Algorithm 910: A Portable C++ Multiple-Precision System for Special-Function Calculations",
+// in ACM TOMS, {VOL 37, ISSUE 4, (February 2011)} (C) ACM, 2011. http://doi.acm.org/10.1145/1916461.1916469
+
+#ifndef _UTIL_FUNCTION_OPERATION_2009_10_27_HPP_
+  #define _UTIL_FUNCTION_OPERATION_2009_10_27_HPP_
+
+  #include <boost/multiprecision/utility/util_function_base.hpp>
+
+  namespace boost
+  {
+    namespace multiprecision
+    {
+      namespace utility
+      {
+        template<typename T> class function_operation : public function_base<T>
+        {
+        protected:
+
+          mutable bool op_ok;
+
+        protected:
+
+          function_operation() : op_ok(false) { }
+
+        public:
+
+          virtual ~function_operation() { }
+
+          bool success(void) const { return op_ok; }
+          T operation(void) const { return my_operation(); }
+
+        protected:
+
+          virtual T my_operation(void) const = 0;
+        };
+      }
+    }
+  }
+
+#endif // _UTIL_FUNCTION_OPERATION_2009_10_27_HPP_

Added: sandbox/multiprecision/boost/utility/util_interpolate.hpp
==============================================================================
--- (empty file)
+++ sandbox/multiprecision/boost/utility/util_interpolate.hpp	2011-09-21 17:00:15 EDT (Wed, 21 Sep 2011)
@@ -0,0 +1,62 @@
+
+//          Copyright Christopher Kormanyos 2002 - 2011.
+// Distributed under the Boost Software License, Version 1.0.
+//    (See accompanying file LICENSE_1_0.txt or copy at
+//          http://www.boost.org/LICENSE_1_0.txt)
+
+// This work is based on an earlier work:
+// "Algorithm 910: A Portable C++ Multiple-Precision System for Special-Function Calculations",
+// in ACM TOMS, {VOL 37, ISSUE 4, (February 2011)} (C) ACM, 2011. http://doi.acm.org/10.1145/1916461.1916469
+
+#ifndef _UTIL_INTERPOLATE_2009_10_27_HPP_
+  #define _UTIL_INTERPOLATE_2009_10_27_HPP_
+
+  #include <vector>
+  #include <algorithm>
+
+  #include <boost/multiprecision/utility/util_point.hpp>
+
+  namespace boost
+  {
+    namespace multiprecision
+    {
+      namespace utility
+      {
+        template<typename T1, typename T2 = T1> struct linear_interpolate
+        {
+          static T2 interpolate(const T1& x, const std::vector<boost::multiprecision::utility::point<T1, T2> >& points)
+          {
+            if(points.empty())
+            {
+              return T2(0);
+            }
+            else if(x <= points.front().x || points.size() == static_cast<std::size_t>(1u))
+            {
+              return points.front().y;
+            }
+            else if(x >= points.back().x)
+            {
+              return points.back().y;
+            }
+            else
+            {
+              const boost::multiprecision::utility::point<T1, T2> x_find(x);
+
+              const typename std::vector<boost::multiprecision::utility::point<T1, T2> >::const_iterator it =
+                std::lower_bound(points.begin(), points.end(), x_find);
+
+              const T1 xn            = (it - 1u)->x;
+              const T1 xnp1_minus_xn = it->x - xn;
+              const T1 delta_x       = x - xn;
+              const T2 yn            = (it - 1u)->y;
+              const T2 ynp1_minus_yn = it->y - yn;
+
+              return T2(yn + T2((delta_x * ynp1_minus_yn) / xnp1_minus_xn));
+            }
+          }
+        };
+      }
+    }
+  }
+
+#endif // _UTIL_INTERPOLATE_2009_10_27_HPP_

Added: sandbox/multiprecision/boost/utility/util_numeric_cast.hpp
==============================================================================
--- (empty file)
+++ sandbox/multiprecision/boost/utility/util_numeric_cast.hpp	2011-09-21 17:00:15 EDT (Wed, 21 Sep 2011)
@@ -0,0 +1,44 @@
+
+//          Copyright Christopher Kormanyos 2002 - 2011.
+// Distributed under the Boost Software License, Version 1.0.
+//    (See accompanying file LICENSE_1_0.txt or copy at
+//          http://www.boost.org/LICENSE_1_0.txt)
+
+// This work is based on an earlier work:
+// "Algorithm 910: A Portable C++ Multiple-Precision System for Special-Function Calculations",
+// in ACM TOMS, {VOL 37, ISSUE 4, (February 2011)} (C) ACM, 2011. http://doi.acm.org/10.1145/1916461.1916469
+
+#ifndef _UTIL_NUMERIC_CAST_2009_11_24_HPP_
+  #define _UTIL_NUMERIC_CAST_2009_11_24_HPP_
+
+  #include <string>
+  #include <sstream>
+
+  namespace boost
+  {
+    namespace multiprecision
+    {
+      namespace utility
+      {
+        template<typename T> inline T numeric_cast(const std::string& str)
+        {
+          std::stringstream ss;
+          ss << str;
+          T t;
+          ss >> t;
+          return t;
+        }
+
+        template<typename T> inline T numeric_cast(const char* const s)
+        {
+          std::stringstream ss;
+          ss << s;
+          T t;
+          ss >> t;
+          return t;
+        }
+      }
+    }
+  }
+
+#endif // _UTIL_NUMERIC_CAST_2009_11_24_HPP_

Added: sandbox/multiprecision/boost/utility/util_point.hpp
==============================================================================
--- (empty file)
+++ sandbox/multiprecision/boost/utility/util_point.hpp	2011-09-21 17:00:15 EDT (Wed, 21 Sep 2011)
@@ -0,0 +1,36 @@
+
+//          Copyright Christopher Kormanyos 2002 - 2011.
+// Distributed under the Boost Software License, Version 1.0.
+//    (See accompanying file LICENSE_1_0.txt or copy at
+//          http://www.boost.org/LICENSE_1_0.txt)
+
+// This work is based on an earlier work:
+// "Algorithm 910: A Portable C++ Multiple-Precision System for Special-Function Calculations",
+// in ACM TOMS, {VOL 37, ISSUE 4, (February 2011)} (C) ACM, 2011. http://doi.acm.org/10.1145/1916461.1916469
+
+#ifndef _UTIL_POINT_2009_10_27_HPP_
+  #define _UTIL_POINT_2009_10_27_HPP_
+
+  namespace boost
+  {
+    namespace multiprecision
+    {
+      namespace utility
+      {
+        template<typename T1, typename T2 = T1> struct point
+        {
+          T1 x;
+          T2 y;
+
+          point(const T1& X = T1(), const T2& Y = T2()) : x(X), y(Y) { }
+        };
+
+        template<typename T1, typename T2> bool inline operator<(const point<T1, T2>& left, const point<T1, T2>& right)
+        {
+          return left.x < right.x;
+        }
+      }
+    }
+  }
+
+#endif // _UTIL_POINT_2009_10_27_HPP_

Added: sandbox/multiprecision/boost/utility/util_power_j_pow_x.hpp
==============================================================================
--- (empty file)
+++ sandbox/multiprecision/boost/utility/util_power_j_pow_x.hpp	2011-09-21 17:00:15 EDT (Wed, 21 Sep 2011)
@@ -0,0 +1,30 @@
+
+//          Copyright Christopher Kormanyos 2002 - 2011.
+// Distributed under the Boost Software License, Version 1.0.
+//    (See accompanying file LICENSE_1_0.txt or copy at
+//          http://www.boost.org/LICENSE_1_0.txt)
+
+// This work is based on an earlier work:
+// "Algorithm 910: A Portable C++ Multiple-Precision System for Special-Function Calculations",
+// in ACM TOMS, {VOL 37, ISSUE 4, (February 2011)} (C) ACM, 2011. http://doi.acm.org/10.1145/1916461.1916469
+
+#ifndef _UTIL_POWER_J_POW_X_2009_01_25_HPP_
+  #define _UTIL_POWER_J_POW_X_2009_01_25_HPP_
+
+  #include <map>
+
+  #include <boost/multiprecision/mp_float_functions.hpp>
+
+  namespace boost
+  {
+    namespace multiprecision
+    {
+      namespace utility
+      {
+        boost::multiprecision::mp_float   j_pow_x(const boost::uint32_t j, const boost::multiprecision::mp_float&   x, std::map<boost::uint32_t, boost::multiprecision::mp_float>&   n_pow_x_prime_factor_map);
+        boost::multiprecision::mp_complex j_pow_x(const boost::uint32_t j, const boost::multiprecision::mp_complex& x, std::map<boost::uint32_t, boost::multiprecision::mp_complex>& n_pow_x_prime_factor_map);
+      }
+    }
+  }
+
+#endif // _UTIL_POWER_J_POW_X_2009_01_25_HPP_

Added: sandbox/multiprecision/boost/utility/util_power_x_pow_n.hpp
==============================================================================
--- (empty file)
+++ sandbox/multiprecision/boost/utility/util_power_x_pow_n.hpp	2011-09-21 17:00:15 EDT (Wed, 21 Sep 2011)
@@ -0,0 +1,72 @@
+
+//          Copyright Christopher Kormanyos 2002 - 2011.
+// Distributed under the Boost Software License, Version 1.0.
+//    (See accompanying file LICENSE_1_0.txt or copy at
+//          http://www.boost.org/LICENSE_1_0.txt)
+
+// This work is based on an earlier work:
+// "Algorithm 910: A Portable C++ Multiple-Precision System for Special-Function Calculations",
+// in ACM TOMS, {VOL 37, ISSUE 4, (February 2011)} (C) ACM, 2011. http://doi.acm.org/10.1145/1916461.1916469
+
+#ifndef _UTIL_POWER_X_POW_N_2009_11_23_HPP_
+  #define _UTIL_POWER_X_POW_N_2009_11_23_HPP_
+
+  namespace boost
+  {
+    namespace multiprecision
+    {
+      namespace utility
+      {
+        template<typename T> inline T x_pow_n_template(const T& t, const boost::int64_t p)
+        {
+          // Compute the pure power of typename T t^p. Binary splitting of the power is
+          // used. The resulting computational complexity has the order of log2[abs(p)].
+
+          if(p < static_cast<boost::int64_t>(0))
+          {
+            return T(1) / x_pow_n_template(t, -p);
+          }
+
+          switch(p)
+          {
+            case static_cast<boost::int64_t>(0):
+              return T(1);
+
+            case static_cast<boost::int64_t>(1):
+              return t;
+
+            case static_cast<boost::int64_t>(2):
+              return t * t;
+
+            case static_cast<boost::int64_t>(3):
+              return (t * t) * t;
+
+            case static_cast<boost::int64_t>(4):
+            {
+              const T tsq = t * t;
+              return tsq * tsq;
+            }
+
+            default:
+            {
+              T value(t);
+              boost::int64_t n;
+
+              for(n = static_cast<boost::int64_t>(1); n <= static_cast<boost::int64_t>(p / static_cast<boost::int64_t>(2)); n *= static_cast<boost::int64_t>(2))
+              {
+                value *= value;
+              }
+
+              const boost::int64_t p_minus_n = static_cast<boost::int64_t>(p - n);
+
+              // Call the function recursively for computing the remaining power of n.
+              return ((p_minus_n == static_cast<boost::int64_t>(0)) ? value
+                                                           : value * boost::multiprecision::utility::x_pow_n_template(t, p_minus_n));
+            }
+          }
+        }
+      }
+    }
+  }
+
+#endif // _UTIL_POWER_X_POW_N_2009_11_23_HPP_

Added: sandbox/multiprecision/boost/utility/util_ranged_function_operation.hpp
==============================================================================
--- (empty file)
+++ sandbox/multiprecision/boost/utility/util_ranged_function_operation.hpp	2011-09-21 17:00:15 EDT (Wed, 21 Sep 2011)
@@ -0,0 +1,46 @@
+
+//          Copyright Christopher Kormanyos 2002 - 2011.
+// Distributed under the Boost Software License, Version 1.0.
+//    (See accompanying file LICENSE_1_0.txt or copy at
+//          http://www.boost.org/LICENSE_1_0.txt)
+
+// This work is based on an earlier work:
+// "Algorithm 910: A Portable C++ Multiple-Precision System for Special-Function Calculations",
+// in ACM TOMS, {VOL 37, ISSUE 4, (February 2011)} (C) ACM, 2011. http://doi.acm.org/10.1145/1916461.1916469
+
+#ifndef _UTIL_RANGED_FUNCTION_OPERATION_2009_10_27_HPP_
+  #define _UTIL_RANGED_FUNCTION_OPERATION_2009_10_27_HPP_
+
+  #include "util_function_operation.hpp"
+
+  namespace boost
+  {
+    namespace multiprecision
+    {
+      namespace utility
+      {
+        template<typename T> class ranged_function_operation : public function_operation<T>
+        {
+        protected:
+
+          const T xlo;
+          const T xhi;
+          const T eps;
+
+        protected:
+
+          ranged_function_operation(const T& lo,
+                                  const T& hi,
+                                  const T& tol) : xlo(lo),
+                                                  xhi(hi),
+                                                  eps(tol) { }
+
+        public:
+
+          virtual ~ranged_function_operation() { }
+        };
+      }
+    }
+  }
+
+#endif // _UTIL_RANGED_FUNCTION_OPERATION_2009_10_27_HPP_

Added: sandbox/multiprecision/boost/utility/util_timer.hpp
==============================================================================
--- (empty file)
+++ sandbox/multiprecision/boost/utility/util_timer.hpp	2011-09-21 17:00:15 EDT (Wed, 21 Sep 2011)
@@ -0,0 +1,45 @@
+
+//          Copyright Christopher Kormanyos 2002 - 2011.
+// Distributed under the Boost Software License, Version 1.0.
+//    (See accompanying file LICENSE_1_0.txt or copy at
+//          http://www.boost.org/LICENSE_1_0.txt)
+
+// This work is based on an earlier work:
+// "Algorithm 910: A Portable C++ Multiple-Precision System for Special-Function Calculations",
+// in ACM TOMS, {VOL 37, ISSUE 4, (February 2011)} (C) ACM, 2011. http://doi.acm.org/10.1145/1916461.1916469
+
+#ifndef _UTIL_TIMER_2010_01_26_HPP_
+  #define _UTIL_TIMER_2010_01_26_HPP_
+
+  #include <boost/noncopyable.hpp>
+
+  namespace boost
+  {
+    namespace multiprecision
+    {
+      namespace utility
+      {
+        struct timer : private boost::noncopyable
+        {
+        public:
+          timer(const double ofs = 0.0) : offset(ofs),
+                                          start (get_sec()) { }
+
+          ~timer() { }
+
+          double elapsed(void) const;
+
+        private:
+          timer(const timer&);
+          const timer& operator=(const timer&);
+
+          const volatile double offset;
+          const volatile double start;
+
+          static double get_sec(void);
+        };
+      }
+    }
+  }
+
+#endif // _UTIL_TIMER_2010_01_26_HPP_

Added: sandbox/multiprecision/boost/utility/util_trapezoid.hpp
==============================================================================
--- (empty file)
+++ sandbox/multiprecision/boost/utility/util_trapezoid.hpp	2011-09-21 17:00:15 EDT (Wed, 21 Sep 2011)
@@ -0,0 +1,81 @@
+
+//          Copyright Christopher Kormanyos 2002 - 2011.
+// Distributed under the Boost Software License, Version 1.0.
+//    (See accompanying file LICENSE_1_0.txt or copy at
+//          http://www.boost.org/LICENSE_1_0.txt)
+
+// This work is based on an earlier work:
+// "Algorithm 910: A Portable C++ Multiple-Precision System for Special-Function Calculations",
+// in ACM TOMS, {VOL 37, ISSUE 4, (February 2011)} (C) ACM, 2011. http://doi.acm.org/10.1145/1916461.1916469
+
+#ifndef _UTIL_TRAPEZOID_2008_09_06_HPP_
+  #define _UTIL_TRAPEZOID_2008_09_06_HPP_
+
+  #include "util_ranged_function_operation.hpp"
+
+  namespace boost
+  {
+    namespace multiprecision
+    {
+      namespace utility
+      {
+        template<typename T> class recursive_trapezoid_rule : public ranged_function_operation<T>
+        {
+        protected:
+
+          recursive_trapezoid_rule(const T& lo, const T& hi, const T& tol) : ranged_function_operation<T>(lo, hi, tol) { }
+
+        public:
+
+          virtual ~recursive_trapezoid_rule() { }
+
+        private:
+
+          virtual T my_operation(void) const
+          {
+            boost::int32_t n = static_cast<boost::int32_t>(1);
+
+            T a = ranged_function_operation<T>::xlo;
+            T b = ranged_function_operation<T>::xhi;
+            T h = (b - a) / T(n);
+        
+            static const T one  = T(1);
+            static const T half = T(0.5);
+
+            T I = (function_base<T>::function(a) + function_base<T>::function(b)) * (h * half);
+
+            for(boost::int32_t k = static_cast<boost::int32_t>(0); k < static_cast<boost::int32_t>(31); k++)
+            {
+              h *= half;
+
+              const T I0 = I;
+
+              T sum(0);
+
+              for(boost::int32_t j = static_cast<boost::int32_t>(1); j <= n; j++)
+              {
+                sum += function(a + (T((j * 2) - 1) * h));
+              }
+
+              I = (I0 * half) + (h * sum);
+
+              const T ratio = I0 / I;
+              const T delta = ((ratio > one) ? (ratio - one) : (one - ratio));
+
+              if((k > static_cast<boost::int32_t>(2)) && delta < ranged_function_operation<T>::eps)
+              {
+                function_operation<T>::op_ok = true;
+                break;
+              }
+
+              n = n * 2;
+            }
+
+            return I;
+          }
+        };
+      }
+    }
+  }
+
+#endif // _UTIL_TRAPEZOID_2008_09_06_HPP_