Import fparserc++_4.5.2.orig.tar.gz

[dgit import orig fparserc++_4.5.2.orig.tar.gz]

4 files changed, 2040 insertions, 0 deletions

diff --git a/mpfr/GmpInt.cc b/mpfr/GmpInt.cc
new file mode 100644
index 0000000..490add4
--- /dev/null
+++ b/mpfr/GmpInt.cc
@@ -0,0 +1,710 @@
+#include "GmpInt.hh"
+#include <gmp.h>
+#include <deque>
+#include <vector>
+#include <cstring>
+#include <cctype>
+
+//===========================================================================
+// Shared data
+//===========================================================================
+namespace
+{
+    unsigned long gIntDefaultNumberOfBits = 256;
+
+    std::vector<char>& intString()
+    {
+        static std::vector<char> str;
+        return str;
+    }
+}
+
+//===========================================================================
+// Auxiliary structs
+//===========================================================================
+struct GmpInt::GmpIntData
+{
+    unsigned mRefCount;
+    GmpIntData* nextFreeNode;
+    mpz_t mInteger;
+
+    GmpIntData(): mRefCount(1), nextFreeNode(0) {}
+};
+
+class GmpInt::GmpIntDataContainer
+{
+    std::deque<GmpInt::GmpIntData> mData;
+    GmpInt::GmpIntData* mFirstFreeNode;
+    GmpInt::GmpIntData* mConst_0;
+
+ public:
+    GmpIntDataContainer(): mFirstFreeNode(0), mConst_0(0) {}
+
+    ~GmpIntDataContainer()
+    {
+        for(size_t i = 0; i < mData.size(); ++i)
+            mpz_clear(mData[i].mInteger);
+    }
+
+    GmpInt::GmpIntData* allocateGmpIntData(unsigned long numberOfBits,
+                                           bool initToZero)
+    {
+        if(mFirstFreeNode)
+        {
+            GmpInt::GmpIntData* node = mFirstFreeNode;
+            mFirstFreeNode = node->nextFreeNode;
+            if(initToZero) mpz_set_si(node->mInteger, 0);
+            ++(node->mRefCount);
+            return node;
+        }
+
+        mData.push_back(GmpInt::GmpIntData());
+        if(numberOfBits > 0)
+            mpz_init2(mData.back().mInteger, numberOfBits);
+        else
+            mpz_init(mData.back().mInteger);
+        return &mData.back();
+    }
+
+    void releaseGmpIntData(GmpIntData* data)
+    {
+        if(--(data->mRefCount) == 0)
+        {
+            data->nextFreeNode = mFirstFreeNode;
+            mFirstFreeNode = data;
+        }
+    }
+
+    GmpInt::GmpIntData* const_0()
+    {
+        if(!mConst_0)
+            mConst_0 = allocateGmpIntData(gIntDefaultNumberOfBits, true);
+        return mConst_0;
+    }
+};
+
+
+GmpInt::GmpIntDataContainer& GmpInt::gmpIntDataContainer()
+{
+    static GmpIntDataContainer container;
+    return container;
+}
+
+//===========================================================================
+// Auxiliary functions
+//===========================================================================
+void GmpInt::setDefaultNumberOfBits(unsigned long value)
+{
+    gIntDefaultNumberOfBits = value;
+}
+
+unsigned long GmpInt::getDefaultNumberOfBits()
+{
+    return gIntDefaultNumberOfBits;
+}
+
+inline void GmpInt::copyIfShared()
+{
+    if(mData->mRefCount > 1)
+    {
+        --(mData->mRefCount);
+        GmpIntData* oldData = mData;
+        mData = gmpIntDataContainer().allocateGmpIntData(0, false);
+        mpz_set(mData->mInteger, oldData->mInteger);
+    }
+}
+
+
+//===========================================================================
+// Constructors, destructor, assignment
+//===========================================================================
+GmpInt::GmpInt(DummyType):
+    mData(gmpIntDataContainer().allocateGmpIntData(0, false))
+{}
+
+GmpInt::GmpInt()
+{
+    mData = gmpIntDataContainer().const_0();
+    ++(mData->mRefCount);
+}
+
+GmpInt::GmpInt(long value)
+{
+    if(value == 0)
+    {
+        mData = gmpIntDataContainer().const_0();
+        ++(mData->mRefCount);
+    }
+    else
+    {
+        mData = gmpIntDataContainer().allocateGmpIntData
+            (gIntDefaultNumberOfBits, false);
+        mpz_set_si(mData->mInteger, value);
+    }
+}
+
+GmpInt::GmpInt(unsigned long value)
+{
+    if(value == 0)
+    {
+        mData = gmpIntDataContainer().const_0();
+        ++(mData->mRefCount);
+    }
+    else
+    {
+        mData = gmpIntDataContainer().allocateGmpIntData
+            (gIntDefaultNumberOfBits, false);
+        mpz_set_ui(mData->mInteger, value);
+    }
+}
+
+GmpInt::GmpInt(int value)
+{
+    if(value == 0)
+    {
+        mData = gmpIntDataContainer().const_0();
+        ++(mData->mRefCount);
+    }
+    else
+    {
+        mData = gmpIntDataContainer().allocateGmpIntData
+            (gIntDefaultNumberOfBits, false);
+        mpz_set_si(mData->mInteger, value);
+    }
+}
+
+GmpInt::GmpInt(double value)
+{
+    const double absValue = value >= 0.0 ? value : -value;
+    if(absValue < 1.0)
+    {
+        mData = gmpIntDataContainer().const_0();
+        ++(mData->mRefCount);
+    }
+    else
+    {
+        mData = gmpIntDataContainer().allocateGmpIntData
+            (gIntDefaultNumberOfBits, false);
+        mpz_set_d(mData->mInteger, value);
+    }
+}
+
+GmpInt::GmpInt(long double value)
+{
+    const long double absValue = value >= 0.0L ? value : -value;
+    if(absValue < 1.0L)
+    {
+        mData = gmpIntDataContainer().const_0();
+        ++(mData->mRefCount);
+    }
+    else
+    {
+        mData = gmpIntDataContainer().allocateGmpIntData
+            (gIntDefaultNumberOfBits, false);
+        mpz_set_d(mData->mInteger, double(value));
+    }
+}
+
+GmpInt::GmpInt(const GmpInt& rhs):
+    mData(rhs.mData)
+{
+    ++(mData->mRefCount);
+}
+
+GmpInt& GmpInt::operator=(const GmpInt& rhs)
+{
+    if(mData != rhs.mData)
+    {
+        gmpIntDataContainer().releaseGmpIntData(mData);
+        mData = rhs.mData;
+        ++(mData->mRefCount);
+    }
+    return *this;
+}
+
+GmpInt& GmpInt::operator=(signed long value)
+{
+    if(value == 0)
+    {
+        gmpIntDataContainer().releaseGmpIntData(mData);
+        mData = gmpIntDataContainer().const_0();
+        ++(mData->mRefCount);
+    }
+    else
+    {
+        if(mData->mRefCount > 1)
+        {
+            --(mData->mRefCount);
+            mData = gmpIntDataContainer().allocateGmpIntData
+                (gIntDefaultNumberOfBits, false);
+        }
+        mpz_set_si(mData->mInteger, value);
+    }
+    return *this;
+}
+
+GmpInt::~GmpInt()
+{
+    gmpIntDataContainer().releaseGmpIntData(mData);
+}
+
+
+//===========================================================================
+// Data getters
+//===========================================================================
+template<>
+void GmpInt::get_raw_mpfr_data<mpz_t>(mpz_t& dest_mpz_t)
+{
+    std::memcpy(&dest_mpz_t, mData->mInteger, sizeof(mpz_t));
+}
+
+const char* GmpInt::getAsString(int base) const
+{
+    intString().resize(mpz_sizeinbase(mData->mInteger, base) + 2);
+    return mpz_get_str(&intString()[0], base, mData->mInteger);
+}
+
+long GmpInt::toInt() const
+{
+    return mpz_get_si(mData->mInteger);
+}
+
+
+//===========================================================================
+// Modifying operators
+//===========================================================================
+GmpInt& GmpInt::operator+=(const GmpInt& rhs)
+{
+    copyIfShared();
+    mpz_add(mData->mInteger, mData->mInteger, rhs.mData->mInteger);
+    return *this;
+}
+
+GmpInt& GmpInt::operator+=(long value)
+{
+    copyIfShared();
+    if(value >= 0)
+        mpz_add_ui(mData->mInteger, mData->mInteger, value);
+    else
+        mpz_sub_ui(mData->mInteger, mData->mInteger, -value);
+    return *this;
+}
+
+GmpInt& GmpInt::operator-=(const GmpInt& rhs)
+{
+    copyIfShared();
+    mpz_sub(mData->mInteger, mData->mInteger, rhs.mData->mInteger);
+    return *this;
+}
+
+GmpInt& GmpInt::operator-=(long value)
+{
+    copyIfShared();
+    if(value >= 0)
+        mpz_sub_ui(mData->mInteger, mData->mInteger, value);
+    else
+        mpz_add_ui(mData->mInteger, mData->mInteger, -value);
+    return *this;
+}
+
+GmpInt& GmpInt::operator*=(const GmpInt& rhs)
+{
+    copyIfShared();
+    mpz_mul(mData->mInteger, mData->mInteger, rhs.mData->mInteger);
+    return *this;
+}
+
+GmpInt& GmpInt::operator*=(long value)
+{
+    copyIfShared();
+    mpz_mul_si(mData->mInteger, mData->mInteger, value);
+    return *this;
+}
+
+GmpInt& GmpInt::operator/=(const GmpInt& rhs)
+{
+    copyIfShared();
+    mpz_tdiv_q(mData->mInteger, mData->mInteger, rhs.mData->mInteger);
+    return *this;
+}
+
+GmpInt& GmpInt::operator/=(long value)
+{
+    copyIfShared();
+    if(value >= 0)
+        mpz_tdiv_q_ui(mData->mInteger, mData->mInteger, value);
+    else
+    {
+        mpz_neg(mData->mInteger, mData->mInteger);
+        mpz_tdiv_q_ui(mData->mInteger, mData->mInteger, -value);
+    }
+    return *this;
+}
+
+GmpInt& GmpInt::operator%=(const GmpInt& rhs)
+{
+    copyIfShared();
+    if(operator<(0))
+    {
+        negate();
+        mpz_mod(mData->mInteger, mData->mInteger, rhs.mData->mInteger);
+        negate();
+    }
+    else
+    {
+        mpz_mod(mData->mInteger, mData->mInteger, rhs.mData->mInteger);
+    }
+    return *this;
+}
+
+GmpInt& GmpInt::operator%=(long value)
+{
+    copyIfShared();
+    if(value < 0) value = -value;
+    if(operator<(0))
+    {
+        negate();
+        mpz_mod_ui(mData->mInteger, mData->mInteger, value);
+        negate();
+    }
+    else
+    {
+        mpz_mod_ui(mData->mInteger, mData->mInteger, value);
+    }
+    return *this;
+}
+
+GmpInt& GmpInt::operator<<=(unsigned long bits)
+{
+    copyIfShared();
+    mpz_mul_2exp(mData->mInteger, mData->mInteger, bits);
+    return *this;
+}
+
+GmpInt& GmpInt::operator>>=(unsigned long bits)
+{
+    copyIfShared();
+    mpz_tdiv_q_2exp(mData->mInteger, mData->mInteger, bits);
+    return *this;
+}
+
+
+//===========================================================================
+// Modifying functions
+//===========================================================================
+void GmpInt::addProduct(const GmpInt& value1, const GmpInt& value2)
+{
+    copyIfShared();
+    mpz_addmul(mData->mInteger, value1.mData->mInteger, value2.mData->mInteger);
+}
+
+void GmpInt::addProduct(const GmpInt& value1, unsigned long value2)
+{
+    copyIfShared();
+    mpz_addmul_ui(mData->mInteger, value1.mData->mInteger, value2);
+}
+
+void GmpInt::subProduct(const GmpInt& value1, const GmpInt& value2)
+{
+    copyIfShared();
+    mpz_submul(mData->mInteger, value1.mData->mInteger, value2.mData->mInteger);
+}
+
+void GmpInt::subProduct(const GmpInt& value1, unsigned long value2)
+{
+    copyIfShared();
+    mpz_submul_ui(mData->mInteger, value1.mData->mInteger, value2);
+}
+
+void GmpInt::negate()
+{
+    copyIfShared();
+    mpz_neg(mData->mInteger, mData->mInteger);
+}
+
+void GmpInt::abs()
+{
+    copyIfShared();
+    mpz_abs(mData->mInteger, mData->mInteger);
+}
+
+GmpInt GmpInt::abs(const GmpInt& value)
+{
+    GmpInt retval(kNoInitialization);
+    mpz_abs(retval.mData->mInteger, value.mData->mInteger);
+    return retval;
+}
+
+
+//===========================================================================
+// Non-modifying operators
+//===========================================================================
+GmpInt GmpInt::operator+(const GmpInt& rhs) const
+{
+    GmpInt retval(kNoInitialization);
+    mpz_add(retval.mData->mInteger, mData->mInteger, rhs.mData->mInteger);
+    return retval;
+}
+
+GmpInt GmpInt::operator+(long value) const
+{
+    GmpInt retval(kNoInitialization);
+    if(value >= 0)
+        mpz_add_ui(retval.mData->mInteger, mData->mInteger, value);
+    else
+        mpz_sub_ui(retval.mData->mInteger, mData->mInteger, -value);
+    return retval;
+}
+
+GmpInt GmpInt::operator-(const GmpInt& rhs) const
+{
+    GmpInt retval(kNoInitialization);
+    mpz_sub(retval.mData->mInteger, mData->mInteger, rhs.mData->mInteger);
+    return retval;
+}
+
+GmpInt GmpInt::operator-(long value) const
+{
+    GmpInt retval(kNoInitialization);
+    if(value >= 0)
+        mpz_sub_ui(retval.mData->mInteger, mData->mInteger, value);
+    else
+        mpz_add_ui(retval.mData->mInteger, mData->mInteger, -value);
+    return retval;
+}
+
+GmpInt GmpInt::operator*(const GmpInt& rhs) const
+{
+    GmpInt retval(kNoInitialization);
+    mpz_mul(retval.mData->mInteger, mData->mInteger, rhs.mData->mInteger);
+    return retval;
+}
+
+GmpInt GmpInt::operator*(long value) const
+{
+    GmpInt retval(kNoInitialization);
+    mpz_mul_si(retval.mData->mInteger, mData->mInteger, value);
+    return retval;
+}
+
+GmpInt GmpInt::operator/(const GmpInt& rhs) const
+{
+    GmpInt retval(kNoInitialization);
+    mpz_tdiv_q(retval.mData->mInteger, mData->mInteger, rhs.mData->mInteger);
+    return retval;
+}
+
+GmpInt GmpInt::operator/(long value) const
+{
+    GmpInt retval(kNoInitialization);
+    if(value >= 0)
+        mpz_tdiv_q_ui(retval.mData->mInteger, mData->mInteger, value);
+    else
+    {
+        mpz_neg(retval.mData->mInteger, mData->mInteger);
+        mpz_tdiv_q_ui(retval.mData->mInteger, retval.mData->mInteger, -value);
+    }
+    return retval;
+}
+
+GmpInt GmpInt::operator%(const GmpInt& rhs) const
+{
+    GmpInt retval(kNoInitialization);
+    if(operator<(0))
+    {
+        mpz_neg(retval.mData->mInteger, mData->mInteger);
+        mpz_mod(retval.mData->mInteger,
+                retval.mData->mInteger, rhs.mData->mInteger);
+        retval.negate();
+    }
+    else
+    {
+        mpz_mod(retval.mData->mInteger, mData->mInteger, rhs.mData->mInteger);
+    }
+    return retval;
+}
+
+GmpInt GmpInt::operator%(long value) const
+{
+    GmpInt retval(kNoInitialization);
+    if(value < 0) value = -value;
+    if(operator<(0))
+    {
+        mpz_neg(retval.mData->mInteger, mData->mInteger);
+        mpz_mod_ui(retval.mData->mInteger, retval.mData->mInteger, value);
+        retval.negate();
+    }
+    else
+    {
+        mpz_mod_ui(retval.mData->mInteger, mData->mInteger, value);
+    }
+    return retval;
+}
+
+GmpInt GmpInt::operator-() const
+{
+    GmpInt retval(kNoInitialization);
+    mpz_neg(retval.mData->mInteger, mData->mInteger);
+    return retval;
+}
+
+GmpInt GmpInt::operator<<(unsigned long bits) const
+{
+    GmpInt retval(kNoInitialization);
+    mpz_mul_2exp(retval.mData->mInteger, mData->mInteger, bits);
+    return retval;
+}
+
+GmpInt GmpInt::operator>>(unsigned long bits) const
+{
+    GmpInt retval(kNoInitialization);
+    mpz_tdiv_q_2exp(retval.mData->mInteger, mData->mInteger, bits);
+    return retval;
+}
+
+
+//===========================================================================
+// Comparison operators
+//===========================================================================
+bool GmpInt::operator<(const GmpInt& rhs) const
+{
+    return mpz_cmp(mData->mInteger, rhs.mData->mInteger) < 0;
+}
+
+bool GmpInt::operator<(long value) const
+{
+    return mpz_cmp_si(mData->mInteger, value) < 0;
+}
+
+bool GmpInt::operator<=(const GmpInt& rhs) const
+{
+    return mpz_cmp(mData->mInteger, rhs.mData->mInteger) <= 0;
+}
+
+bool GmpInt::operator<=(long value) const
+{
+    return mpz_cmp_si(mData->mInteger, value) <= 0;
+}
+
+bool GmpInt::operator>(const GmpInt& rhs) const
+{
+    return mpz_cmp(mData->mInteger, rhs.mData->mInteger) > 0;
+}
+
+bool GmpInt::operator>(long value) const
+{
+    return mpz_cmp_si(mData->mInteger, value) > 0;
+}
+
+bool GmpInt::operator>=(const GmpInt& rhs) const
+{
+    return mpz_cmp(mData->mInteger, rhs.mData->mInteger) >= 0;
+}
+
+bool GmpInt::operator>=(long value) const
+{
+    return mpz_cmp_si(mData->mInteger, value) >= 0;
+}
+
+bool GmpInt::operator==(const GmpInt& rhs) const
+{
+    return mpz_cmp(mData->mInteger, rhs.mData->mInteger) == 0;
+}
+
+bool GmpInt::operator==(long value) const
+{
+    return mpz_cmp_si(mData->mInteger, value) == 0;
+}
+
+bool GmpInt::operator!=(const GmpInt& rhs) const
+{
+    return mpz_cmp(mData->mInteger, rhs.mData->mInteger) != 0;
+}
+
+bool GmpInt::operator!=(long value) const
+{
+    return mpz_cmp_si(mData->mInteger, value) != 0;
+}
+
+void GmpInt::parseValue(const char* value)
+{
+    mpz_set_str(mData->mInteger, value, 10);
+}
+
+void GmpInt::parseValue(const char* value, char** endptr)
+{
+    static std::vector<char> str;
+
+    unsigned startIndex = 0;
+    while(value[startIndex] && std::isspace(value[startIndex])) ++startIndex;
+    if(!value[startIndex]) { *endptr = const_cast<char*>(value); return; }
+
+    unsigned endIndex = startIndex;
+    if(value[endIndex] == '-') ++endIndex;
+    if(!std::isdigit(value[endIndex]))
+    { *endptr = const_cast<char*>(value); return; }
+    if(value[endIndex] == '0' && value[endIndex+1] == 'x')
+    {
+        endIndex += 1;
+        while(std::isxdigit(value[++endIndex])) {}
+    }
+    else
+    {
+        while(std::isdigit(value[++endIndex])) {}
+    }
+
+    str.reserve(endIndex - startIndex + 1);
+    str.assign(value + startIndex, value + endIndex);
+    str.push_back(0);
+
+    mpz_set_str(mData->mInteger, &str[0], 0);
+    *endptr = const_cast<char*>(value + endIndex);
+}
+
+GmpInt GmpInt::parseString(const char* str, char** endptr)
+{
+    GmpInt retval(kNoInitialization);
+    retval.parseValue(str, endptr);
+    return retval;
+}
+
+//===========================================================================
+// Operator functions
+//===========================================================================
+GmpInt operator+(long lhs, const GmpInt& rhs)
+{
+    GmpInt retval(GmpInt::kNoInitialization);
+    if(lhs >= 0)
+        mpz_add_ui(retval.mData->mInteger, rhs.mData->mInteger, lhs);
+    else
+        mpz_sub_ui(retval.mData->mInteger, rhs.mData->mInteger, -lhs);
+    return retval;
+}
+
+GmpInt operator-(long lhs, const GmpInt& rhs)
+{
+    GmpInt retval(GmpInt::kNoInitialization);
+    if(lhs >= 0)
+        mpz_ui_sub(retval.mData->mInteger, lhs, rhs.mData->mInteger);
+    else
+    {
+        mpz_add_ui(retval.mData->mInteger, rhs.mData->mInteger, -lhs);
+        mpz_neg(retval.mData->mInteger, retval.mData->mInteger);
+    }
+    return retval;
+}
+
+GmpInt operator*(long lhs, const GmpInt& rhs)
+{
+    return rhs * lhs;
+}
+
+GmpInt operator/(long lhs, const GmpInt& rhs)
+{
+    return GmpInt(lhs) / rhs;
+}
+
+GmpInt operator%(long lhs, const GmpInt& rhs)
+{
+    return GmpInt(lhs) % rhs;
+}
diff --git a/mpfr/GmpInt.hh b/mpfr/GmpInt.hh
new file mode 100644
index 0000000..1c1c171
--- /dev/null
+++ b/mpfr/GmpInt.hh
@@ -0,0 +1,148 @@
+#ifndef ONCE_FP_GMP_INT_HH_
+#define ONCE_FP_GMP_INT_HH_
+
+#include <iostream>
+
+class GmpInt
+{
+ public:
+    /* A default of 256 bits will be used for all newly-instantiated GmpInt
+       objects. This default can be changed with the function below.
+    */
+    static void setDefaultNumberOfBits(unsigned long);
+    static unsigned long getDefaultNumberOfBits();
+
+    GmpInt();
+    GmpInt(long value);
+    GmpInt(unsigned long value);
+    GmpInt(int value);
+    GmpInt(double value);
+    GmpInt(long double value);
+
+    GmpInt(const GmpInt&);
+    GmpInt& operator=(const GmpInt&);
+    GmpInt& operator=(signed long value);
+
+    ~GmpInt();
+
+
+    /* This function can be used to retrieve the raw mpz_t data structure
+       used by this object. (The template trick is used to avoid a dependency
+       of this header file with <gmp.h>.)
+       In other words, it can be called like:
+
+         mpz_t raw_mpz_data;
+         intValue.get_raw_mpz_data(raw_mpz_data);
+
+       Note that the returned mpz_t should be considered as read-only and
+       not be modified from the outside because it may be shared among
+       several objects. If the calling code needs to modify the data, it
+       should copy it for itself first with the appropriate GMP library
+       functions.
+     */
+    template<typename Mpz_t>
+    void get_raw_mpfr_data(Mpz_t& dest_mpz_t);
+
+
+    // Note that the returned char* points to an internal (shared) buffer
+    // which will be valid until the next time this function is called
+    // (by any object).
+    const char* getAsString(int base = 10) const;
+    long toInt() const;
+
+    GmpInt& operator+=(const GmpInt&);
+    GmpInt& operator+=(long);
+    GmpInt& operator-=(const GmpInt&);
+    GmpInt& operator-=(long);
+    GmpInt& operator*=(const GmpInt&);
+    GmpInt& operator*=(long);
+    GmpInt& operator/=(const GmpInt&);
+    GmpInt& operator/=(long);
+    GmpInt& operator%=(const GmpInt&);
+    GmpInt& operator%=(long);
+
+    GmpInt& operator<<=(unsigned long);
+    GmpInt& operator>>=(unsigned long);
+
+    // Equivalent to "+= value1 * value2;"
+    void addProduct(const GmpInt& value1, const GmpInt& value2);
+    void addProduct(const GmpInt& value1, unsigned long value2);
+
+    // Equivalent to "-= value1 * value2;"
+    void subProduct(const GmpInt& value1, const GmpInt& value2);
+    void subProduct(const GmpInt& value1, unsigned long value2);
+
+    void negate();
+    void abs();
+    static GmpInt abs(const GmpInt&);
+
+    GmpInt operator+(const GmpInt&) const;
+    GmpInt operator+(long) const;
+    GmpInt operator-(const GmpInt&) const;
+    GmpInt operator-(long) const;
+    GmpInt operator*(const GmpInt&) const;
+    GmpInt operator*(long) const;
+    GmpInt operator/(const GmpInt&) const;
+    GmpInt operator/(long) const;
+    GmpInt operator%(const GmpInt&) const;
+    GmpInt operator%(long) const;
+
+    GmpInt operator-() const;
+
+    GmpInt operator<<(unsigned long) const;
+    GmpInt operator>>(unsigned long) const;
+
+    bool operator<(const GmpInt&) const;
+    bool operator<(long) const;
+    bool operator<=(const GmpInt&) const;
+    bool operator<=(long) const;
+    bool operator>(const GmpInt&) const;
+    bool operator>(long) const;
+    bool operator>=(const GmpInt&) const;
+    bool operator>=(long) const;
+    bool operator==(const GmpInt&) const;
+    bool operator==(long) const;
+    bool operator!=(const GmpInt&) const;
+    bool operator!=(long) const;
+
+    void parseValue(const char* value);
+    void parseValue(const char* value, char** endptr);
+    static GmpInt parseString(const char* str, char** endptr);
+
+
+ private:
+    struct GmpIntData;
+    class GmpIntDataContainer;
+
+    GmpIntData* mData;
+
+    enum DummyType { kNoInitialization };
+    GmpInt(DummyType);
+
+    void copyIfShared();
+    static GmpIntDataContainer& gmpIntDataContainer();
+
+    friend GmpInt operator+(long lhs, const GmpInt& rhs);
+    friend GmpInt operator-(long lhs, const GmpInt& rhs);
+};
+
+GmpInt operator+(long lhs, const GmpInt& rhs);
+GmpInt operator-(long lhs, const GmpInt& rhs);
+GmpInt operator*(long lhs, const GmpInt& rhs);
+GmpInt operator/(long lhs, const GmpInt& rhs);
+GmpInt operator%(long lhs, const GmpInt& rhs);
+
+inline bool operator<(long lhs, const GmpInt& rhs) { return rhs > lhs; }
+inline bool operator<=(long lhs, const GmpInt& rhs) { return rhs >= lhs; }
+inline bool operator>(long lhs, const GmpInt& rhs) { return rhs < lhs; }
+inline bool operator>=(long lhs, const GmpInt& rhs) { return rhs <= lhs; }
+inline bool operator==(long lhs, const GmpInt& rhs) { return rhs == lhs; }
+inline bool operator!=(long lhs, const GmpInt& rhs) { return rhs != lhs; }
+
+inline std::ostream& operator<<(std::ostream& os, const GmpInt& value)
+{
+    os << value.getAsString();
+    return os;
+}
+
+#endif
diff --git a/mpfr/MpfrFloat.cc b/mpfr/MpfrFloat.cc
new file mode 100644
index 0000000..112c684
--- /dev/null
+++ b/mpfr/MpfrFloat.cc
@@ -0,0 +1,976 @@
+#include "MpfrFloat.hh"
+#include <stdio.h>
+#include <mpfr.h>
+#include <deque>
+#include <vector>
+#include <cstring>
+#include <cassert>
+
+//===========================================================================
+// Auxiliary structs
+//===========================================================================
+struct MpfrFloat::MpfrFloatData
+{
+    unsigned mRefCount;
+    MpfrFloatData* nextFreeNode;
+    mpfr_t mFloat;
+
+    MpfrFloatData(): mRefCount(1), nextFreeNode(0) {}
+};
+
+class MpfrFloat::MpfrFloatDataContainer
+{
+    unsigned long mDefaultPrecision;
+    std::deque<MpfrFloatData> mData;
+    MpfrFloatData* mFirstFreeNode;
+
+    MpfrFloatData
+    *mConst_0, *mConst_pi, *mConst_e, *mConst_log2, *mConst_epsilon;
+
+    void recalculateEpsilon()
+    {
+        mpfr_set_si(mConst_epsilon->mFloat, 1, GMP_RNDN);
+        mpfr_div_2ui(mConst_epsilon->mFloat, mConst_epsilon->mFloat,
+                     mDefaultPrecision*7/8 - 1, GMP_RNDN);
+    }
+
+ public:
+    MpfrFloatDataContainer():
+        mDefaultPrecision(256), mFirstFreeNode(0), mConst_0(0),
+        mConst_pi(0), mConst_e(0), mConst_log2(0), mConst_epsilon(0)
+    {}
+
+    ~MpfrFloatDataContainer()
+    {
+        for(size_t i = 0; i < mData.size(); ++i)
+            mpfr_clear(mData[i].mFloat);
+    }
+
+    MpfrFloatData* allocateMpfrFloatData(bool initToZero)
+    {
+        if(mFirstFreeNode)
+        {
+            MpfrFloatData* node = mFirstFreeNode;
+            mFirstFreeNode = node->nextFreeNode;
+            if(initToZero) mpfr_set_si(node->mFloat, 0, GMP_RNDN);
+            ++(node->mRefCount);
+            return node;
+        }
+
+        mData.push_back(MpfrFloatData());
+        mpfr_init2(mData.back().mFloat, mDefaultPrecision);
+        if(initToZero) mpfr_set_si(mData.back().mFloat, 0, GMP_RNDN);
+        return &mData.back();
+    }
+
+    void releaseMpfrFloatData(MpfrFloatData* data)
+    {
+        if(--(data->mRefCount) == 0)
+        {
+            data->nextFreeNode = mFirstFreeNode;
+            mFirstFreeNode = data;
+        }
+    }
+
+    void setDefaultPrecision(unsigned long bits)
+    {
+        if(bits != mDefaultPrecision)
+        {
+            mDefaultPrecision = bits;
+            for(size_t i = 0; i < mData.size(); ++i)
+                mpfr_prec_round(mData[i].mFloat, bits, GMP_RNDN);
+
+            if(mConst_pi) mpfr_const_pi(mConst_pi->mFloat, GMP_RNDN);
+            if(mConst_e)
+            {
+                mpfr_set_si(mConst_e->mFloat, 1, GMP_RNDN);
+                mpfr_exp(mConst_e->mFloat, mConst_e->mFloat, GMP_RNDN);
+            }
+            if(mConst_log2) mpfr_const_log2(mConst_log2->mFloat, GMP_RNDN);
+            if(mConst_epsilon) recalculateEpsilon();
+        }
+    }
+
+    unsigned long getDefaultPrecision() const
+    {
+        return mDefaultPrecision;
+    }
+
+    MpfrFloatData* const_0()
+    {
+        if(!mConst_0) mConst_0 = allocateMpfrFloatData(true);
+        return mConst_0;
+    }
+
+    MpfrFloat const_pi()
+    {
+        if(!mConst_pi)
+        {
+            mConst_pi = allocateMpfrFloatData(false);
+            mpfr_const_pi(mConst_pi->mFloat, GMP_RNDN);
+        }
+        return MpfrFloat(mConst_pi);
+    }
+
+    MpfrFloat const_e()
+    {
+        if(!mConst_e)
+        {
+            mConst_e = allocateMpfrFloatData(false);
+            mpfr_set_si(mConst_e->mFloat, 1, GMP_RNDN);
+            mpfr_exp(mConst_e->mFloat, mConst_e->mFloat, GMP_RNDN);
+        }
+        return MpfrFloat(mConst_e);
+    }
+
+    MpfrFloat const_log2()
+    {
+        if(!mConst_log2)
+        {
+            mConst_log2 = allocateMpfrFloatData(false);
+            mpfr_const_log2(mConst_log2->mFloat, GMP_RNDN);
+        }
+        return MpfrFloat(mConst_log2);
+    }
+
+    MpfrFloat const_epsilon()
+    {
+        if(!mConst_epsilon)
+        {
+            mConst_epsilon = allocateMpfrFloatData(false);
+            recalculateEpsilon();
+        }
+        return MpfrFloat(mConst_epsilon);
+    }
+};
+
+
+//===========================================================================
+// Shared data
+//===========================================================================
+// This should ensure that the container is not accessed by any MpfrFloat
+// instance before it has been constructed or after it has been destroyed
+// (which might otherwise happen if MpfrFloat is instantiated globally.)
+MpfrFloat::MpfrFloatDataContainer& MpfrFloat::mpfrFloatDataContainer()
+{
+    static MpfrFloat::MpfrFloatDataContainer container;
+    return container;
+}
+
+
+//===========================================================================
+// Auxiliary functions
+//===========================================================================
+void MpfrFloat::setDefaultMantissaBits(unsigned long bits)
+{
+    mpfrFloatDataContainer().setDefaultPrecision(bits);
+}
+
+unsigned long MpfrFloat::getCurrentDefaultMantissaBits()
+{
+    return mpfrFloatDataContainer().getDefaultPrecision();
+}
+
+inline void MpfrFloat::copyIfShared()
+{
+    if(mData->mRefCount > 1)
+    {
+        --(mData->mRefCount);
+        MpfrFloatData* oldData = mData;
+        mData = mpfrFloatDataContainer().allocateMpfrFloatData(false);
+        mpfr_set(mData->mFloat, oldData->mFloat, GMP_RNDN);
+    }
+}
+
+
+//===========================================================================
+// Constructors, destructor, assignment
+//===========================================================================
+MpfrFloat::MpfrFloat(DummyType):
+    mData(mpfrFloatDataContainer().allocateMpfrFloatData(false))
+{}
+
+MpfrFloat::MpfrFloat(MpfrFloatData* data):
+    mData(data)
+{
+    assert(data != 0);
+    ++(mData->mRefCount);
+}
+
+MpfrFloat::MpfrFloat():
+    mData(mpfrFloatDataContainer().const_0())
+{
+    ++(mData->mRefCount);
+}
+
+MpfrFloat::MpfrFloat(double value)
+{
+    if(value == 0.0)
+    {
+        mData = mpfrFloatDataContainer().const_0();
+        ++(mData->mRefCount);
+    }
+    else
+    {
+        mData = mpfrFloatDataContainer().allocateMpfrFloatData(false);
+        mpfr_set_d(mData->mFloat, value, GMP_RNDN);
+    }
+}
+
+MpfrFloat::MpfrFloat(long double value)
+{
+    if(value == 0.0L)
+    {
+        mData = mpfrFloatDataContainer().const_0();
+        ++(mData->mRefCount);
+    }
+    else
+    {
+        mData = mpfrFloatDataContainer().allocateMpfrFloatData(false);
+        mpfr_set_ld(mData->mFloat, value, GMP_RNDN);
+    }
+}
+
+MpfrFloat::MpfrFloat(long value)
+{
+    if(value == 0)
+    {
+        mData = mpfrFloatDataContainer().const_0();
+        ++(mData->mRefCount);
+    }
+    else
+    {
+        mData = mpfrFloatDataContainer().allocateMpfrFloatData(false);
+        mpfr_set_si(mData->mFloat, value, GMP_RNDN);
+    }
+}
+
+MpfrFloat::MpfrFloat(int value)
+{
+    if(value == 0)
+    {
+        mData = mpfrFloatDataContainer().const_0();
+        ++(mData->mRefCount);
+    }
+    else
+    {
+        mData = mpfrFloatDataContainer().allocateMpfrFloatData(false);
+        mpfr_set_si(mData->mFloat, value, GMP_RNDN);
+    }
+}
+
+MpfrFloat::MpfrFloat(const char* value, char** endptr):
+    mData(mpfrFloatDataContainer().allocateMpfrFloatData(false))
+{
+    mpfr_strtofr(mData->mFloat, value, endptr, 0, GMP_RNDN);
+}
+
+MpfrFloat::~MpfrFloat()
+{
+    mpfrFloatDataContainer().releaseMpfrFloatData(mData);
+}
+
+MpfrFloat::MpfrFloat(const MpfrFloat& rhs):
+    mData(rhs.mData)
+{
+    ++(mData->mRefCount);
+}
+
+MpfrFloat& MpfrFloat::operator=(const MpfrFloat& rhs)
+{
+    if(mData != rhs.mData)
+    {
+        mpfrFloatDataContainer().releaseMpfrFloatData(mData);
+        mData = rhs.mData;
+        ++(mData->mRefCount);
+    }
+    return *this;
+}
+
+MpfrFloat& MpfrFloat::operator=(double value)
+{
+    if(value == 0.0)
+    {
+        mpfrFloatDataContainer().releaseMpfrFloatData(mData);
+        mData = mpfrFloatDataContainer().const_0();
+        ++(mData->mRefCount);
+    }
+    else
+    {
+        if(mData->mRefCount > 1)
+        {
+            --(mData->mRefCount);
+            mData = mpfrFloatDataContainer().allocateMpfrFloatData(false);
+        }
+        mpfr_set_d(mData->mFloat, value, GMP_RNDN);
+    }
+    return *this;
+}
+
+MpfrFloat& MpfrFloat::operator=(long double value)
+{
+    if(value == 0.0L)
+    {
+        mpfrFloatDataContainer().releaseMpfrFloatData(mData);
+        mData = mpfrFloatDataContainer().const_0();
+        ++(mData->mRefCount);
+    }
+    else
+    {
+        if(mData->mRefCount > 1)
+        {
+            --(mData->mRefCount);
+            mData = mpfrFloatDataContainer().allocateMpfrFloatData(false);
+        }
+        mpfr_set_ld(mData->mFloat, value, GMP_RNDN);
+    }
+    return *this;
+}
+
+MpfrFloat& MpfrFloat::operator=(long value)
+{
+    if(value == 0)
+    {
+        mpfrFloatDataContainer().releaseMpfrFloatData(mData);
+        mData = mpfrFloatDataContainer().const_0();
+        ++(mData->mRefCount);
+    }
+    else
+    {
+        if(mData->mRefCount > 1)
+        {
+            --(mData->mRefCount);
+            mData = mpfrFloatDataContainer().allocateMpfrFloatData(false);
+        }
+        mpfr_set_si(mData->mFloat, value, GMP_RNDN);
+    }
+    return *this;
+}
+
+MpfrFloat& MpfrFloat::operator=(int value)
+{
+    if(value == 0)
+    {
+        mpfrFloatDataContainer().releaseMpfrFloatData(mData);
+        mData = mpfrFloatDataContainer().const_0();
+        ++(mData->mRefCount);
+    }
+    else
+    {
+        if(mData->mRefCount > 1)
+        {
+            --(mData->mRefCount);
+            mData = mpfrFloatDataContainer().allocateMpfrFloatData(false);
+        }
+        mpfr_set_si(mData->mFloat, value, GMP_RNDN);
+    }
+    return *this;
+}
+
+/*
+MpfrFloat& MpfrFloat::operator=(const char* value)
+{
+    if(mData->mRefCount > 1)
+    {
+        --(mData->mRefCount);
+        mData = mpfrFloatDataContainer().allocateMpfrFloatData(false);
+    }
+
+    mpfr_set_str(mData->mFloat, value, 10, GMP_RNDN);
+    return *this;
+}
+*/
+
+void MpfrFloat::parseValue(const char* value)
+{
+    copyIfShared();
+    mpfr_set_str(mData->mFloat, value, 10, GMP_RNDN);
+}
+
+void MpfrFloat::parseValue(const char* value, char** endptr)
+{
+    copyIfShared();
+    mpfr_strtofr(mData->mFloat, value, endptr, 0, GMP_RNDN);
+}
+
+
+//===========================================================================
+// Data getters
+//===========================================================================
+template<>
+void MpfrFloat::get_raw_mpfr_data<mpfr_t>(mpfr_t& dest_mpfr_t)
+{
+    std::memcpy(&dest_mpfr_t, mData->mFloat, sizeof(mpfr_t));
+}
+
+const char* MpfrFloat::getAsString(unsigned precision) const
+{
+#if(MPFR_VERSION_MAJOR < 2 || (MPFR_VERSION_MAJOR == 2 && MPFR_VERSION_MINOR < 4))
+    static const char* retval =
+        "[mpfr_snprintf() is not supported in mpfr versions prior to 2.4]";
+    return retval;
+#else
+    static std::vector<char> str;
+    str.resize(precision+30);
+    mpfr_snprintf(&(str[0]), precision+30, "%.*RNg", precision, mData->mFloat);
+    return &(str[0]);
+#endif
+}
+
+bool MpfrFloat::isInteger() const
+{
+    return mpfr_integer_p(mData->mFloat) != 0;
+}
+
+long MpfrFloat::toInt() const
+{
+    return mpfr_get_si(mData->mFloat, GMP_RNDN);
+}
+
+double MpfrFloat::toDouble() const
+{
+    return mpfr_get_d(mData->mFloat, GMP_RNDN);
+}
+
+
+//===========================================================================
+// Modifying operators
+//===========================================================================
+MpfrFloat& MpfrFloat::operator+=(const MpfrFloat& rhs)
+{
+    copyIfShared();
+    mpfr_add(mData->mFloat, mData->mFloat, rhs.mData->mFloat, GMP_RNDN);
+    return *this;
+}
+
+MpfrFloat& MpfrFloat::operator+=(double value)
+{
+    copyIfShared();
+    mpfr_add_d(mData->mFloat, mData->mFloat, value, GMP_RNDN);
+    return *this;
+}
+
+MpfrFloat& MpfrFloat::operator-=(const MpfrFloat& rhs)
+{
+    copyIfShared();
+    mpfr_sub(mData->mFloat, mData->mFloat, rhs.mData->mFloat, GMP_RNDN);
+    return *this;
+}
+
+MpfrFloat& MpfrFloat::operator-=(double value)
+{
+    copyIfShared();
+    mpfr_sub_d(mData->mFloat, mData->mFloat, value, GMP_RNDN);
+    return *this;
+}
+
+MpfrFloat& MpfrFloat::operator*=(const MpfrFloat& rhs)
+{
+    copyIfShared();
+    mpfr_mul(mData->mFloat, mData->mFloat, rhs.mData->mFloat, GMP_RNDN);
+    return *this;
+}
+
+MpfrFloat& MpfrFloat::operator*=(double value)
+{
+    copyIfShared();
+    mpfr_mul_d(mData->mFloat, mData->mFloat, value, GMP_RNDN);
+    return *this;
+}
+
+MpfrFloat& MpfrFloat::operator/=(const MpfrFloat& rhs)
+{
+    copyIfShared();
+    mpfr_div(mData->mFloat, mData->mFloat, rhs.mData->mFloat, GMP_RNDN);
+    return *this;
+}
+
+MpfrFloat& MpfrFloat::operator/=(double value)
+{
+    copyIfShared();
+    mpfr_div_d(mData->mFloat, mData->mFloat, value, GMP_RNDN);
+    return *this;
+}
+
+MpfrFloat& MpfrFloat::operator%=(const MpfrFloat& rhs)
+{
+    copyIfShared();
+    mpfr_fmod(mData->mFloat, mData->mFloat, rhs.mData->mFloat, GMP_RNDN);
+    return *this;
+}
+
+
+//===========================================================================
+// Modifying functions
+//===========================================================================
+void MpfrFloat::negate()
+{
+    copyIfShared();
+    mpfr_neg(mData->mFloat, mData->mFloat, GMP_RNDN);
+}
+
+void MpfrFloat::abs()
+{
+    copyIfShared();
+    mpfr_abs(mData->mFloat, mData->mFloat, GMP_RNDN);
+}
+
+
+//===========================================================================
+// Non-modifying operators
+//===========================================================================
+MpfrFloat MpfrFloat::operator+(const MpfrFloat& rhs) const
+{
+    MpfrFloat retval(kNoInitialization);
+    mpfr_add(retval.mData->mFloat, mData->mFloat, rhs.mData->mFloat, GMP_RNDN);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::operator+(double value) const
+{
+    MpfrFloat retval(kNoInitialization);
+    mpfr_add_d(retval.mData->mFloat, mData->mFloat, value, GMP_RNDN);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::operator-(const MpfrFloat& rhs) const
+{
+    MpfrFloat retval(kNoInitialization);
+    mpfr_sub(retval.mData->mFloat, mData->mFloat, rhs.mData->mFloat, GMP_RNDN);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::operator-(double value) const
+{
+    MpfrFloat retval(kNoInitialization);
+    mpfr_sub_d(retval.mData->mFloat, mData->mFloat, value, GMP_RNDN);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::operator*(const MpfrFloat& rhs) const
+{
+    MpfrFloat retval(kNoInitialization);
+    mpfr_mul(retval.mData->mFloat, mData->mFloat, rhs.mData->mFloat, GMP_RNDN);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::operator*(double value) const
+{
+    MpfrFloat retval(kNoInitialization);
+    mpfr_mul_d(retval.mData->mFloat, mData->mFloat, value, GMP_RNDN);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::operator/(const MpfrFloat& rhs) const
+{
+    MpfrFloat retval(kNoInitialization);
+    mpfr_div(retval.mData->mFloat, mData->mFloat, rhs.mData->mFloat, GMP_RNDN);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::operator/(double value) const
+{
+    MpfrFloat retval(kNoInitialization);
+    mpfr_div_d(retval.mData->mFloat, mData->mFloat, value, GMP_RNDN);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::operator%(const MpfrFloat& rhs) const
+{
+    MpfrFloat retval(kNoInitialization);
+    mpfr_fmod(retval.mData->mFloat, mData->mFloat, rhs.mData->mFloat, GMP_RNDN);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::operator-() const
+{
+    MpfrFloat retval(kNoInitialization);
+    mpfr_neg(retval.mData->mFloat, mData->mFloat, GMP_RNDN);
+    return retval;
+}
+
+
+
+//===========================================================================
+// Comparison operators
+//===========================================================================
+bool MpfrFloat::operator<(const MpfrFloat& rhs) const
+{
+    return mpfr_cmp(mData->mFloat, rhs.mData->mFloat) < 0;
+}
+
+bool MpfrFloat::operator<(double value) const
+{
+    return mpfr_cmp_d(mData->mFloat, value) < 0;
+}
+
+bool MpfrFloat::operator<=(const MpfrFloat& rhs) const
+{
+    return mpfr_cmp(mData->mFloat, rhs.mData->mFloat) <= 0;
+}
+
+bool MpfrFloat::operator<=(double value) const
+{
+    return mpfr_cmp_d(mData->mFloat, value) <= 0;
+}
+
+bool MpfrFloat::operator>(const MpfrFloat& rhs) const
+{
+    return mpfr_cmp(mData->mFloat, rhs.mData->mFloat) > 0;
+}
+
+bool MpfrFloat::operator>(double value) const
+{
+    return mpfr_cmp_d(mData->mFloat, value) > 0;
+}
+
+bool MpfrFloat::operator>=(const MpfrFloat& rhs) const
+{
+    return mpfr_cmp(mData->mFloat, rhs.mData->mFloat) >= 0;
+}
+
+bool MpfrFloat::operator>=(double value) const
+{
+    return mpfr_cmp_d(mData->mFloat, value) >= 0;
+}
+
+bool MpfrFloat::operator==(const MpfrFloat& rhs) const
+{
+    return mpfr_cmp(mData->mFloat, rhs.mData->mFloat) == 0;
+}
+
+bool MpfrFloat::operator==(double value) const
+{
+    return mpfr_cmp_d(mData->mFloat, value) == 0;
+}
+
+bool MpfrFloat::operator!=(const MpfrFloat& rhs) const
+{
+    return mpfr_cmp(mData->mFloat, rhs.mData->mFloat) != 0;
+}
+
+bool MpfrFloat::operator!=(double value) const
+{
+    return mpfr_cmp_d(mData->mFloat, value) != 0;
+}
+
+
+//===========================================================================
+// Operator functions
+//===========================================================================
+MpfrFloat operator+(double lhs, const MpfrFloat& rhs)
+{
+    MpfrFloat retval(MpfrFloat::kNoInitialization);
+    mpfr_add_d(retval.mData->mFloat, rhs.mData->mFloat, lhs, GMP_RNDN);
+    return retval;
+}
+
+MpfrFloat operator-(double lhs, const MpfrFloat& rhs)
+{
+    MpfrFloat retval(MpfrFloat::kNoInitialization);
+    mpfr_d_sub(retval.mData->mFloat, lhs, rhs.mData->mFloat, GMP_RNDN);
+    return retval;
+}
+
+MpfrFloat operator*(double lhs, const MpfrFloat& rhs)
+{
+    return rhs * lhs;
+}
+
+MpfrFloat operator/(double lhs, const MpfrFloat& rhs)
+{
+    return MpfrFloat(lhs) / rhs;
+}
+
+MpfrFloat operator%(double lhs, const MpfrFloat& rhs)
+{
+    return MpfrFloat(lhs) % rhs;
+}
+
+std::ostream& operator<<(std::ostream& os, const MpfrFloat& value)
+{
+    os << value.getAsString(unsigned(os.precision()));
+    return os;
+}
+
+//===========================================================================
+// Static functions
+//===========================================================================
+MpfrFloat MpfrFloat::log(const MpfrFloat& value)
+{
+    MpfrFloat retval(MpfrFloat::kNoInitialization);
+    mpfr_log(retval.mData->mFloat, value.mData->mFloat, GMP_RNDN);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::log2(const MpfrFloat& value)
+{
+    MpfrFloat retval(MpfrFloat::kNoInitialization);
+    mpfr_log2(retval.mData->mFloat, value.mData->mFloat, GMP_RNDN);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::log10(const MpfrFloat& value)
+{
+    MpfrFloat retval(MpfrFloat::kNoInitialization);
+    mpfr_log10(retval.mData->mFloat, value.mData->mFloat, GMP_RNDN);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::exp(const MpfrFloat& value)
+{
+    MpfrFloat retval(MpfrFloat::kNoInitialization);
+    mpfr_exp(retval.mData->mFloat, value.mData->mFloat, GMP_RNDN);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::exp2(const MpfrFloat& value)
+{
+    MpfrFloat retval(MpfrFloat::kNoInitialization);
+    mpfr_exp2(retval.mData->mFloat, value.mData->mFloat, GMP_RNDN);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::exp10(const MpfrFloat& value)
+{
+    MpfrFloat retval(MpfrFloat::kNoInitialization);
+    mpfr_exp10(retval.mData->mFloat, value.mData->mFloat, GMP_RNDN);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::cos(const MpfrFloat& value)
+{
+    MpfrFloat retval(MpfrFloat::kNoInitialization);
+    mpfr_cos(retval.mData->mFloat, value.mData->mFloat, GMP_RNDN);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::sin(const MpfrFloat& value)
+{
+    MpfrFloat retval(MpfrFloat::kNoInitialization);
+    mpfr_sin(retval.mData->mFloat, value.mData->mFloat, GMP_RNDN);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::tan(const MpfrFloat& value)
+{
+    MpfrFloat retval(MpfrFloat::kNoInitialization);
+    mpfr_tan(retval.mData->mFloat, value.mData->mFloat, GMP_RNDN);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::sec(const MpfrFloat& value)
+{
+    MpfrFloat retval(MpfrFloat::kNoInitialization);
+    mpfr_sec(retval.mData->mFloat, value.mData->mFloat, GMP_RNDN);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::csc(const MpfrFloat& value)
+{
+    MpfrFloat retval(MpfrFloat::kNoInitialization);
+    mpfr_csc(retval.mData->mFloat, value.mData->mFloat, GMP_RNDN);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::cot(const MpfrFloat& value)
+{
+    MpfrFloat retval(MpfrFloat::kNoInitialization);
+    mpfr_cot(retval.mData->mFloat, value.mData->mFloat, GMP_RNDN);
+    return retval;
+}
+
+void MpfrFloat::sincos(const MpfrFloat& value,
+                       MpfrFloat& sin,
+                       MpfrFloat& cos)
+{
+    sin.copyIfShared();
+    cos.copyIfShared();
+    mpfr_sin_cos
+        (sin.mData->mFloat, cos.mData->mFloat, value.mData->mFloat, GMP_RNDN);
+}
+
+MpfrFloat MpfrFloat::acos(const MpfrFloat& value)
+{
+    MpfrFloat retval(MpfrFloat::kNoInitialization);
+    mpfr_acos(retval.mData->mFloat, value.mData->mFloat, GMP_RNDN);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::asin(const MpfrFloat& value)
+{
+    MpfrFloat retval(MpfrFloat::kNoInitialization);
+    mpfr_asin(retval.mData->mFloat, value.mData->mFloat, GMP_RNDN);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::atan(const MpfrFloat& value)
+{
+    MpfrFloat retval(MpfrFloat::kNoInitialization);
+    mpfr_atan(retval.mData->mFloat, value.mData->mFloat, GMP_RNDN);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::atan2(const MpfrFloat& value1, const MpfrFloat& value2)
+{
+    MpfrFloat retval(MpfrFloat::kNoInitialization);
+    mpfr_atan2(retval.mData->mFloat,
+               value1.mData->mFloat, value2.mData->mFloat, GMP_RNDN);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::hypot(const MpfrFloat& value1, const MpfrFloat& value2)
+{
+    MpfrFloat retval(MpfrFloat::kNoInitialization);
+    mpfr_hypot(retval.mData->mFloat,
+               value1.mData->mFloat, value2.mData->mFloat, GMP_RNDN);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::cosh(const MpfrFloat& value)
+{
+    MpfrFloat retval(MpfrFloat::kNoInitialization);
+    mpfr_cosh(retval.mData->mFloat, value.mData->mFloat, GMP_RNDN);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::sinh(const MpfrFloat& value)
+{
+    MpfrFloat retval(MpfrFloat::kNoInitialization);
+    mpfr_sinh(retval.mData->mFloat, value.mData->mFloat, GMP_RNDN);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::tanh(const MpfrFloat& value)
+{
+    MpfrFloat retval(MpfrFloat::kNoInitialization);
+    mpfr_tanh(retval.mData->mFloat, value.mData->mFloat, GMP_RNDN);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::acosh(const MpfrFloat& value)
+{
+    MpfrFloat retval(MpfrFloat::kNoInitialization);
+    mpfr_acosh(retval.mData->mFloat, value.mData->mFloat, GMP_RNDN);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::asinh(const MpfrFloat& value)
+{
+    MpfrFloat retval(MpfrFloat::kNoInitialization);
+    mpfr_asinh(retval.mData->mFloat, value.mData->mFloat, GMP_RNDN);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::atanh(const MpfrFloat& value)
+{
+    MpfrFloat retval(MpfrFloat::kNoInitialization);
+    mpfr_atanh(retval.mData->mFloat, value.mData->mFloat, GMP_RNDN);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::sqrt(const MpfrFloat& value)
+{
+    MpfrFloat retval(MpfrFloat::kNoInitialization);
+    mpfr_sqrt(retval.mData->mFloat, value.mData->mFloat, GMP_RNDN);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::cbrt(const MpfrFloat& value)
+{
+    MpfrFloat retval(MpfrFloat::kNoInitialization);
+    mpfr_cbrt(retval.mData->mFloat, value.mData->mFloat, GMP_RNDN);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::root(const MpfrFloat& value, unsigned long root)
+{
+    MpfrFloat retval(MpfrFloat::kNoInitialization);
+    mpfr_root(retval.mData->mFloat, value.mData->mFloat, root, GMP_RNDN);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::pow(const MpfrFloat& value1, const MpfrFloat& value2)
+{
+    MpfrFloat retval(MpfrFloat::kNoInitialization);
+    mpfr_pow(retval.mData->mFloat,
+             value1.mData->mFloat, value2.mData->mFloat, GMP_RNDN);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::pow(const MpfrFloat& value, long exponent)
+{
+    MpfrFloat retval(MpfrFloat::kNoInitialization);
+    mpfr_pow_si(retval.mData->mFloat, value.mData->mFloat, exponent, GMP_RNDN);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::abs(const MpfrFloat& value)
+{
+    MpfrFloat retval(MpfrFloat::kNoInitialization);
+    mpfr_abs(retval.mData->mFloat, value.mData->mFloat, GMP_RNDN);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::dim(const MpfrFloat& value1, const MpfrFloat& value2)
+{
+    MpfrFloat retval(MpfrFloat::kNoInitialization);
+    mpfr_dim(retval.mData->mFloat,
+             value1.mData->mFloat, value2.mData->mFloat, GMP_RNDN);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::round(const MpfrFloat& value)
+{
+    MpfrFloat retval(MpfrFloat::kNoInitialization);
+    mpfr_round(retval.mData->mFloat, value.mData->mFloat);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::ceil(const MpfrFloat& value)
+{
+    MpfrFloat retval(MpfrFloat::kNoInitialization);
+    mpfr_ceil(retval.mData->mFloat, value.mData->mFloat);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::floor(const MpfrFloat& value)
+{
+    MpfrFloat retval(MpfrFloat::kNoInitialization);
+    mpfr_floor(retval.mData->mFloat, value.mData->mFloat);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::trunc(const MpfrFloat& value)
+{
+    MpfrFloat retval(MpfrFloat::kNoInitialization);
+    mpfr_trunc(retval.mData->mFloat, value.mData->mFloat);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::parseString(const char* str, char** endptr)
+{
+    MpfrFloat retval(MpfrFloat::kNoInitialization);
+    mpfr_strtofr(retval.mData->mFloat, str, endptr, 0, GMP_RNDN);
+    return retval;
+}
+
+MpfrFloat MpfrFloat::const_pi()
+{
+    return mpfrFloatDataContainer().const_pi();
+}
+
+MpfrFloat MpfrFloat::const_e()
+{
+    return mpfrFloatDataContainer().const_e();
+}
+
+MpfrFloat MpfrFloat::const_log2()
+{
+    return mpfrFloatDataContainer().const_log2();
+}
+
+MpfrFloat MpfrFloat::someEpsilon()
+{
+    return mpfrFloatDataContainer().const_epsilon();
+}
diff --git a/mpfr/MpfrFloat.hh b/mpfr/MpfrFloat.hh
new file mode 100644
index 0000000..d455d24
--- /dev/null
+++ b/mpfr/MpfrFloat.hh
@@ -0,0 +1,206 @@
+#ifndef ONCE_FP_MPFR_FLOAT_
+#define ONCE_FP_MPFR_FLOAT_
+
+#include <iostream>
+
+class MpfrFloat
+{
+ public:
+    /* A default of 256 bits will be used unless changed with this function.
+       Note that all existing and cached GMP objects will be resized to the
+       specified precision (which can be a somewhat heavy operation).
+    */
+    static void setDefaultMantissaBits(unsigned long bits);
+
+    static unsigned long getCurrentDefaultMantissaBits();
+
+    /* The default constructor initializes the object to the value 0.
+       It's efficient to instantiate such zero-initialized objects because
+       all of them will share the same mpfr data. (Also any object initialized
+       with or assigned the explicit value of zero will also share that one
+       mpfr data.) Thus multiple zero-initialized MpfrFloat instances won't
+       consume significant amounts of memory (until they are modified to
+       contain some other value, of course).
+
+       Important caveat:
+       ----------------
+       Note that initializing an MpfrFloat object with, for example, 0.1 will
+       suffer from accuracy problems (at least if the MpfrFloat object has
+       more mantissa bits than a double). The C++ double value 0.1 has only
+       53 mantissa bits, while the MpfrFloat object usually has more. If the
+       MpfrFloat object is initialized with a double, only that many bits of
+       accuracy will end up in the value of the MpfrFloat object. This can
+       create significant rounding/accuracy problems in some cases.
+       If you need to initialize the MpfrObject with some value (which cannot
+       be represented accurately by base-2 floating point numbers, eg. 0.1)
+       at full mantissa precision, you have to use parseValue("0.1") instead,
+       rather than relying on the constructor taking a double type value.
+    */
+    MpfrFloat();
+    MpfrFloat(double value);
+    MpfrFloat(long double value);
+    MpfrFloat(long value);
+    MpfrFloat(int value);
+    MpfrFloat(const char* value, char** endptr);
+
+    ~MpfrFloat();
+
+    MpfrFloat(const MpfrFloat&);
+
+    MpfrFloat& operator=(const MpfrFloat&);
+    MpfrFloat& operator=(double value);
+    MpfrFloat& operator=(long double value);
+    MpfrFloat& operator=(long value);
+    MpfrFloat& operator=(int value);
+    //MpfrFloat& operator=(const char* value);
+
+    void parseValue(const char* value);
+    void parseValue(const char* value, char** endptr);
+
+
+    /* This function can be used to retrieve the raw mpfr_t data structure
+       used by this object. (The template trick is used to avoid a dependency
+       of this header file with <mpfr.h>.)
+       In other words, it can be called like:
+
+         mpfr_t raw_mpfr_data;
+         floatValue.get_raw_mpfr_data(raw_mpfr_data);
+
+       Note that the returned mpfr_t should be considered as read-only and
+       not be modified from the outside because it may be shared among
+       several objects. If the calling code needs to modify the data, it
+       should copy it for itself first with the appropriate MPFR library
+       functions.
+     */
+    template<typename Mpfr_t>
+    void get_raw_mpfr_data(Mpfr_t& dest_mpfr_t);
+
+
+    /* Note that the returned char* points to an internal (shared) buffer
+       which will be valid until the next time this function is called
+       (by any object).
+    */
+    const char* getAsString(unsigned precision) const;
+
+    bool isInteger() const;
+    long toInt() const;
+    double toDouble() const;
+
+    MpfrFloat& operator+=(const MpfrFloat&);
+    MpfrFloat& operator+=(double);
+    MpfrFloat& operator-=(const MpfrFloat&);
+    MpfrFloat& operator-=(double);
+    MpfrFloat& operator*=(const MpfrFloat&);
+    MpfrFloat& operator*=(double);
+    MpfrFloat& operator/=(const MpfrFloat&);
+    MpfrFloat& operator/=(double);
+    MpfrFloat& operator%=(const MpfrFloat&);
+
+    void negate();
+    void abs();
+
+    MpfrFloat operator+(const MpfrFloat&) const;
+    MpfrFloat operator+(double) const;
+    MpfrFloat operator-(const MpfrFloat&) const;
+    MpfrFloat operator-(double) const;
+    MpfrFloat operator*(const MpfrFloat&) const;
+    MpfrFloat operator*(double) const;
+    MpfrFloat operator/(const MpfrFloat&) const;
+    MpfrFloat operator/(double) const;
+    MpfrFloat operator%(const MpfrFloat&) const;
+
+    MpfrFloat operator-() const;
+
+    bool operator<(const MpfrFloat&) const;
+    bool operator<(double) const;
+    bool operator<=(const MpfrFloat&) const;
+    bool operator<=(double) const;
+    bool operator>(const MpfrFloat&) const;
+    bool operator>(double) const;
+    bool operator>=(const MpfrFloat&) const;
+    bool operator>=(double) const;
+    bool operator==(const MpfrFloat&) const;
+    bool operator==(double) const;
+    bool operator!=(const MpfrFloat&) const;
+    bool operator!=(double) const;
+
+    static MpfrFloat log(const MpfrFloat&);
+    static MpfrFloat log2(const MpfrFloat&);
+    static MpfrFloat log10(const MpfrFloat&);
+    static MpfrFloat exp(const MpfrFloat&);
+    static MpfrFloat exp2(const MpfrFloat&);
+    static MpfrFloat exp10(const MpfrFloat&);
+    static MpfrFloat cos(const MpfrFloat&);
+    static MpfrFloat sin(const MpfrFloat&);
+    static MpfrFloat tan(const MpfrFloat&);
+    static MpfrFloat sec(const MpfrFloat&);
+    static MpfrFloat csc(const MpfrFloat&);
+    static MpfrFloat cot(const MpfrFloat&);
+    static void sincos(const MpfrFloat&, MpfrFloat& sin, MpfrFloat& cos);
+    static MpfrFloat acos(const MpfrFloat&);
+    static MpfrFloat asin(const MpfrFloat&);
+    static MpfrFloat atan(const MpfrFloat&);
+    static MpfrFloat atan2(const MpfrFloat&, const MpfrFloat&);
+    static MpfrFloat hypot(const MpfrFloat&, const MpfrFloat&);
+    static MpfrFloat cosh(const MpfrFloat&);
+    static MpfrFloat sinh(const MpfrFloat&);
+    static MpfrFloat tanh(const MpfrFloat&);
+    static MpfrFloat acosh(const MpfrFloat&);
+    static MpfrFloat asinh(const MpfrFloat&);
+    static MpfrFloat atanh(const MpfrFloat&);
+    static MpfrFloat sqrt(const MpfrFloat&);
+    static MpfrFloat cbrt(const MpfrFloat&);
+    static MpfrFloat root(const MpfrFloat&, unsigned long root);
+    static MpfrFloat pow(const MpfrFloat&, const MpfrFloat&);
+    static MpfrFloat pow(const MpfrFloat&, long exponent);
+    static MpfrFloat abs(const MpfrFloat&);
+    static MpfrFloat dim(const MpfrFloat&, const MpfrFloat&);
+    static MpfrFloat round(const MpfrFloat&);
+    static MpfrFloat ceil(const MpfrFloat&);
+    static MpfrFloat floor(const MpfrFloat&);
+    static MpfrFloat trunc(const MpfrFloat&);
+
+    static MpfrFloat parseString(const char* str, char** endptr);
+
+    // These values are cached (and recalculated every time the mantissa bits
+    // change), so it's efficient to call these repeatedly:
+    static MpfrFloat const_pi();
+    static MpfrFloat const_e();
+    static MpfrFloat const_log2();
+    static MpfrFloat someEpsilon();
+
+
+ private:
+    struct MpfrFloatData;
+    class MpfrFloatDataContainer;
+
+    MpfrFloatData* mData;
+
+    enum DummyType { kNoInitialization };
+    MpfrFloat(DummyType);
+    MpfrFloat(MpfrFloatData*);
+
+    void copyIfShared();
+    static MpfrFloatDataContainer& mpfrFloatDataContainer();
+
+    friend MpfrFloat operator+(double lhs, const MpfrFloat& rhs);
+    friend MpfrFloat operator-(double lhs, const MpfrFloat& rhs);
+};
+
+MpfrFloat operator+(double lhs, const MpfrFloat& rhs);
+MpfrFloat operator-(double lhs, const MpfrFloat& rhs);
+MpfrFloat operator*(double lhs, const MpfrFloat& rhs);
+MpfrFloat operator/(double lhs, const MpfrFloat& rhs);
+MpfrFloat operator%(double lhs, const MpfrFloat& rhs);
+
+inline bool operator<(double lhs, const MpfrFloat& rhs) { return rhs > lhs; }
+inline bool operator<=(double lhs, const MpfrFloat& rhs) { return rhs >= lhs; }
+inline bool operator>(double lhs, const MpfrFloat& rhs) { return rhs < lhs; }
+inline bool operator>=(double lhs, const MpfrFloat& rhs) { return rhs <= lhs; }
+inline bool operator==(double lhs, const MpfrFloat& rhs) { return rhs == lhs; }
+inline bool operator!=(double lhs, const MpfrFloat& rhs) { return rhs != lhs; }
+
+// This function takes into account the value of os.precision()
+std::ostream& operator<<(std::ostream& os, const MpfrFloat& value);
+
+#endif