X-Git-Url: https://git.ucc.asn.au/?p=ipdf%2Fcode.git;a=blobdiff_plain;f=src%2Frational.h;h=41cce09078b4e96eb17a6bb6e28478ea64637b57;hp=b3709ad7f47b067f30986c90ab09f00bb17d8d4f;hb=00ef152dc3a065e37fe45c1cd8023d739f518b8e;hpb=33356addacfe4296ecb613c6c4696f082e351159 diff --git a/src/rational.h b/src/rational.h index b3709ad..41cce09 100644 --- a/src/rational.h +++ b/src/rational.h @@ -12,8 +12,7 @@ namespace IPDF { -/** Greatest Common Divisor - Euclid's algorithm **/ - +/* Recursive version of GCD template T gcd(const T & a, const T & b) { @@ -24,12 +23,12 @@ T gcd(const T & a, const T & b) if (a > b) return gcd(a-b,b); return gcd(a, b-a); } +*/ -/* +/** Greatest Common Divisor of p and q **/ template T gcd(const T & p, const T & q) { - Debug("p/q = % T g(1); T big(p); T small(q); @@ -47,15 +46,15 @@ T gcd(const T & p, const T & q) } return small; } -*/ + template struct Rational { /** Construct from a double.**/ - Rational(double d = 0) : P(d*1e3), Q(1e3) // Possibly the worst thing ever... + Rational(double d=0) : P(d*1e6), Q(1e6) // Possibly the worst thing ever... { Simplify(); - //CheckAccuracy(d, "Construct from double"); + CheckAccuracy(d, "Construct from double"); } Rational(const T & _P, const T & _Q) : P(_P), Q(_Q) @@ -75,7 +74,11 @@ struct Rational P = (P < 0) ? -P : P; Q = -Q; } - + if (P == 0) + { + Q = 1; + return; + } T g = gcd(llabs(P),llabs(Q)); P /= g; Q /= g; @@ -94,9 +97,6 @@ struct Rational bool operator>=(const Rational & r) const {return *this == r || *this > r;} bool operator!=(const Rational & r) const {return !(*this == r);} - - - /* Rational operator+(const Rational & r) const { Rational result = (r.P == 0) ? Rational(P,Q) : Rational(P*r.Q + r.P*Q, Q*r.Q); @@ -109,7 +109,6 @@ struct Rational result.CheckAccuracy(ToDouble() - r.ToDouble(),"-"); return result; } - */ Rational operator*(const Rational & r) const { Rational result(P * r.P, Q * r.Q); @@ -121,16 +120,17 @@ struct Rational } Rational operator/(const Rational & r) const { - Rational result = (r.P == 0) ? Rational(P,Q) : Rational(P*r.Q + r.P*Q, Q*r.Q); - if (!result.CheckAccuracy(ToDouble() / r.ToDouble(),"*")) + 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()); } return result; } - Rational operator+(const Rational & r) const {return Rational(ToDouble()+r.ToDouble());} - Rational operator-(const Rational & r) const {return Rational(ToDouble()-r.ToDouble());} + /** To cheat, use these **/ + //Rational operator+(const Rational & r) const {return Rational(ToDouble()+r.ToDouble());} + //Rational operator-(const Rational & r) const {return Rational(ToDouble()-r.ToDouble());} //Rational operator*(const Rational & r) const {return Rational(ToDouble()*r.ToDouble());} //Rational operator/(const Rational & r) const {return Rational(ToDouble()/r.ToDouble());} @@ -141,7 +141,7 @@ struct Rational Rational & operator/=(const Rational & r) {this->operator=(*this/r); return *this;} double ToDouble() const {return (double)(P) / (double)(Q);} - bool CheckAccuracy(double d, const char * msg, double threshold = 1e-6) const + bool CheckAccuracy(double d, const char * msg, double threshold = 1e-3) const { double result = fabs(ToDouble() - d) / d; if (result > threshold)