src/gmprat.h

   1 /**
   2  * @file gmpint.h
   3  * @brief Wraps to GMP mpq_t type using inlines
   4  */
   5
   6 #ifndef _GMPRAT_H
   7 #define _GMPRAT_H
   8
   9 #include <gmp.h>
  10 #include <string>
  11
  12 class Gmprat
  13 {
  14         public:
  15                 Gmprat(double d=0) {mpq_init(m_op); mpq_set_d(m_op, d);}
  16                 //Gmprat(int64_t p = 0, int64_t q = 1) {mpq_init(m_op); mpq_set_si(m_op, p, q);}
  17                 //Gmprat(const std::string & str, int base=10) {mpq_init_set_str(m_op, str.c_str(), base);}
  18                 Gmprat(const Gmprat & cpy) {mpq_init(m_op); mpq_set(m_op, cpy.m_op);}
  19                 virtual ~Gmprat() {mpq_clear(m_op);}
  20
  21                 //operator int64_t() const {return mpq_get_si(m_op);}
  22                 //operator uint64_t() const {return mpq_get_ui(m_op);}
  23                 operator double() const {return mpq_get_d(m_op);}
  24                 //operator float() const {return (float)ToDouble();}
  25                 double ToDouble() const {return mpq_get_d(m_op);}
  26                 std::string Str() const
  27                 {
  28                         //TODO: Make less hacky, if we care.
  29                         // Convert to scientific notation
  30                         if (Negative())
  31                         {
  32                                 return "-" + operator-().Str();
  33                         }
  34                         double p = Log10();
  35                         if (isinf(p))
  36                                 return "0";
  37
  38                         int P = (int)p;
  39                         double C = pow(10.0, p - P);
  40                         std::stringstream s;
  41                         s << C;
  42                         if (P != 0)
  43                                 s << "e"<< P;
  44                         //s << "("<<ToDouble()<<")";
  45                         return s.str();
  46                 }
  47
  48                 double Log10() const
  49                 {
  50                         mpz_t num; mpz_init(num); mpq_get_num(num, m_op);
  51                         mpz_t den; mpz_init(den); mpq_get_den(den, m_op);
  52
  53                         double lognum = 0;
  54                         double logden = 0;
  55                         while (mpz_sizeinbase(num, 10) > 10)
  56                         {
  57
  58                                 mpz_div_ui(num, num, 1e10);
  59                                 lognum += 10;
  60                         }
  61                         uint64_t n = mpz_get_ui(num);
  62                         if (n == 0)
  63                         {
  64                                 return -INFINITY;
  65                         }
  66                         lognum += log(n)/log(10.0);
  67                         //Debug("%lu", mpz_get_ui(den));
  68                         while (mpz_sizeinbase(den, 10) > 10)
  69                         {
  70                                 mpz_div_ui(den, den, 1e10);
  71                                 logden += 10;
  72                         }
  73                         uint64_t d = mpz_get_ui(den);
  74                         // if d is zero, its been rounded down we hope
  75                         if (d != 0)
  76                                 logden += log(d)/log(10.0);
  77
  78                         return (lognum - logden);
  79                 }
  80
  81
  82                 bool Negative() const {return (mpz_sgn(mpq_numref(m_op)) < 0);}
  83
  84                 Gmprat & operator=(const Gmprat & equ) {mpq_set(m_op, equ.m_op); return *this;}
  85                 Gmprat & operator=(const double & equ) {mpq_set_d(m_op, equ); return *this;}
  86                 Gmprat & operator+=(const Gmprat & add) {mpq_add(m_op, m_op, add.m_op); return *this;}
  87                 Gmprat & operator-=(const Gmprat & sub) {mpq_sub(m_op, m_op, sub.m_op); return *this;}
  88                 Gmprat & operator*=(const Gmprat & mul) {mpq_mul(m_op, m_op, mul.m_op); return *this;}
  89                 Gmprat & operator/=(const Gmprat & div) {mpq_div(m_op, m_op, div.m_op); return *this;}
  90
  91                 Gmprat operator+(const Gmprat & add) const {Gmprat a(*this); a += add; return a;}
  92                 Gmprat operator-(const Gmprat & sub) const {Gmprat a(*this); a -= sub; return a;}
  93                 Gmprat operator*(const Gmprat & mul) const {Gmprat a(*this); a *= mul; return a;}
  94                 Gmprat operator/(const Gmprat & div) const {Gmprat a(*this); a /= div; return a;}
  95                 //Gmprat operator%(const Gmprat & div) const {Gmprat a(*this); mpq_mod(a.m_op, a.m_op, div.m_op); return a;}
  96                 Gmprat operator-() const {return (Gmprat(0L)-*this);}
  97
  98
  99                 bool operator==(const Gmprat & cmp) const {return mpq_cmp(m_op, cmp.m_op) == 0;}
 100                 bool operator!=(const Gmprat & cmp) const {return mpq_cmp(m_op, cmp.m_op) != 0;}
 101                 bool operator<(const Gmprat & cmp) const {return mpq_cmp(m_op, cmp.m_op) < 0;}
 102                 bool operator>(const Gmprat & cmp) const {return mpq_cmp(m_op, cmp.m_op) > 0;}
 103                 bool operator<=(const Gmprat & cmp) const {return mpq_cmp(m_op, cmp.m_op) <= 0;}
 104                 bool operator>=(const Gmprat & cmp) const {return mpq_cmp(m_op, cmp.m_op) >= 0;}
 105
 106                 Gmprat Abs() const {Gmprat a(*this); mpq_abs(a.m_op, a.m_op); return a;}
 107
 108
 109         private:
 110                 friend std::ostream& operator<<(std::ostream& os, const Gmprat & fith);
 111                 mpq_t m_op;
 112 };
 113
 114 inline std::ostream & operator<<(std::ostream & os, const Gmprat & fith)
 115 {
 116         os << fith.Str();
 117         return os;
 118 }
 119
 120
 121
 122
 123 #endif //_GMPRAT_H
 124