X-Git-Url: https://git.ucc.asn.au/?p=ipdf%2Fcode.git;a=blobdiff_plain;f=src%2Frational.h;h=e74567db3529cd2e2ee520eac3e5c4eb70558268;hp=6269ef181c3c6ad7a7126ab17f743a077487448f;hb=a182cbba4fe9e6a36c5063735dbec1c5340da04c;hpb=a7da45c62d5a0604785479f5dfeb1c2a65d0f9de diff --git a/src/rational.h b/src/rational.h index 6269ef1..e74567d 100644 --- a/src/rational.h +++ b/src/rational.h @@ -8,9 +8,17 @@ #include "common.h" #include #include +#include "arbint.h" namespace IPDF { + +template T Tabs(const T & a) +{ + return llabs(a); +} + + /* Recursive version of GCD template @@ -65,7 +73,10 @@ struct Rational Rational(double d=0) : P(d*1e6), Q(1e6) // Possibly the worst thing ever... { Simplify(); - CheckAccuracy(d, "Construct from double"); + //if (!CheckAccuracy(d, "Construct from double")) + { + //Fatal("Bwah bwah :("); + } } Rational(const T & _P, const T & _Q) : P(_P), Q(_Q) @@ -90,8 +101,8 @@ struct Rational Q = T(1); return; } - T g = gcd(T(llabs(P)),T(llabs(Q))); - Debug("Got gcd!"); + T g = gcd(Tabs(P), Tabs(Q)); + //Debug("Got gcd!"); P /= g; Q /= g; } @@ -112,7 +123,10 @@ struct Rational Rational operator+(const Rational & r) const { Rational result = (r.P == T(0)) ? Rational(P,Q) : Rational(P*r.Q + r.P*Q, Q*r.Q); - result.CheckAccuracy(ToDouble() + r.ToDouble(),"+"); + //if (!result.CheckAccuracy(ToDouble() * r.ToDouble(),"+")) + //{ + // Debug("This is %s (%f) and r is %s (%f)", Str().c_str(), ToDouble(), r.Str().c_str(), r.ToDouble()); + //} return result; } Rational operator-(const Rational & r) const @@ -124,19 +138,19 @@ struct Rational Rational operator*(const Rational & r) const { Rational result(P * r.P, Q * r.Q); - if (!result.CheckAccuracy(ToDouble() * r.ToDouble(),"*")) - { - Debug("This is %s (%f) and r is %s (%f)", Str().c_str(), ToDouble(), r.Str().c_str(), r.ToDouble()); - } + //if (!result.CheckAccuracy(ToDouble() * r.ToDouble(),"*")) + //{ + // Debug("This is %s (%f) and r is %s (%f)", Str().c_str(), ToDouble(), r.Str().c_str(), r.ToDouble()); + //} return result; } Rational operator/(const Rational & r) const { Rational result(P * r.Q, Q*r.P); - if (!result.CheckAccuracy(ToDouble() / r.ToDouble(),"/")) - { - Debug("This is %s (%f) and r is %s (%f)", Str().c_str(), ToDouble(), r.Str().c_str(), r.ToDouble()); - } + //if (!result.CheckAccuracy(ToDouble() / r.ToDouble(),"/")) + //{ + // Debug("This is %s (%f) and r is %s (%f)", Str().c_str(), ToDouble(), r.Str().c_str(), r.ToDouble()); + //} return result; } @@ -155,10 +169,12 @@ struct Rational double ToDouble() const {return (double)(P) / (double)(Q);} bool CheckAccuracy(double d, const char * msg, double threshold = 1e-3) const { - double result = fabs(ToDouble() - d) / d; + double result = fabs(ToDouble() - d); + if (d != 0e0) result /= d; if (result > threshold) { Warn("(%s) : Rational %s (%f) is not close enough at representing %f (%f vs %f)", msg, Str().c_str(), ToDouble(), d, result, threshold); + Backtrace(); return false; } return true; @@ -166,7 +182,7 @@ struct Rational std::string Str() const { std::stringstream s; - s << (int64_t)P << "/" << (int64_t)Q; + s << int64_t(P) << "/" << int64_t(Q); return s.str(); } @@ -184,6 +200,7 @@ inline Rational pow(const Rational & a, const Rational