diff options
author | Norbert Preining <norbert@preining.info> | 2019-09-02 13:46:59 +0900 |
---|---|---|
committer | Norbert Preining <norbert@preining.info> | 2019-09-02 13:46:59 +0900 |
commit | e0c6872cf40896c7be36b11dcc744620f10adf1d (patch) | |
tree | 60335e10d2f4354b0674ec22d7b53f0f8abee672 /graphics/epix/frac.cc |
Initial commit
Diffstat (limited to 'graphics/epix/frac.cc')
-rw-r--r-- | graphics/epix/frac.cc | 262 |
1 files changed, 262 insertions, 0 deletions
diff --git a/graphics/epix/frac.cc b/graphics/epix/frac.cc new file mode 100644 index 0000000000..6fadd246af --- /dev/null +++ b/graphics/epix/frac.cc @@ -0,0 +1,262 @@ +/* + * frac.cc -- ePiX rational number class + * + * This file is part of ePiX, a C++ library for creating high-quality + * figures in LaTeX + * + * Version 1.2.0-2 + * Last Change: September 26, 2007 + */ + +/* + * Copyright (C) 2001, 2002, 2003, 2004, 2005, 2006, 2007 + * Andrew D. Hwang <rot 13 nujnat at zngupf dot ubylpebff dot rqh> + * Department of Mathematics and Computer Science + * College of the Holy Cross + * Worcester, MA, 01610-2395, USA + */ + +/* + * ePiX is free software; you can redistribute it and/or modify it + * under the terms of the GNU General Public License as published by + * the Free Software Foundation; either version 2 of the License, or + * (at your option) any later version. + * + * ePiX is distributed in the hope that it will be useful, but WITHOUT + * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY + * or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public + * License for more details. + * + * You should have received a copy of the GNU General Public License + * along with ePiX; if not, write to the Free Software Foundation, Inc., + * 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA + */ +#include <cmath> + +#include "frac.h" + +namespace ePiX { + + const int MAX_DENOM(10000); + const double EPS(1.0/MAX_DENOM); + + unsigned int __epix_gcd (int i, unsigned int j); + + frac::frac(const int n, const unsigned int d) + : m_num(n), m_denom(d) { } + + // express arg as a fraction with denominator no larger than MAX_DENOM + frac::frac(double arg) + { + unsigned int denom(1); + if (fabs(arg) < EPS) + { + m_num = 0; + m_denom = 1; + } + + else + { + int sgn(arg < 0 ? -1 : 1); + double abs_arg(sgn*arg); + + // store best approximation + int best_num(0); + double running_error(0.5); + + while (denom <= 1+MAX_DENOM) + { + const double tmp(abs_arg*denom); + + int tmp_lo((int) floor(tmp)); + + if (tmp - tmp_lo < EPS) // good approx + { + best_num = tmp_lo; + break; + } + + else if (tmp - tmp_lo < running_error) // improved approx + { + best_num = tmp_lo; + running_error = tmp - best_num; + } + + int tmp_hi((int) ceil(tmp)); + if (tmp_hi - tmp < EPS) + { + best_num = tmp_hi; + break; + } + + else if (tmp_hi - tmp < running_error) + { + best_num = tmp_hi; + running_error = best_num - tmp; + } + + ++denom; + } + + m_num = sgn*best_num; + m_denom = denom; + } + } // end of frac::frac(double arg) + + + // increment operators + frac& frac::operator += (const frac& arg) + { + m_num *= arg.m_denom; + m_num += m_denom*arg.m_num; + m_denom *= arg.m_denom; + + return *this; + } + + frac& frac::operator -= (const frac& arg) + { + m_num *= arg.m_denom; + m_num -= m_denom*arg.m_num; + m_denom *= arg.m_denom; + + return *this; + } + + frac& frac::operator *= (const frac& arg) + { + m_num *= arg.m_num; + m_denom *= arg.m_denom; + + return *this; + } + + frac& frac::operator /= (const frac& arg) + { + unsigned int arg_num(arg.m_num < 0 ? -arg.m_num : arg.m_num); + + m_num *= arg.m_denom; + m_denom *= arg_num; + + if (arg.m_num < 0) + m_num = -1; + + return *this; + } + + frac& frac::reduce() + { + unsigned int factor(__epix_gcd(m_num, m_denom)); + m_num /= factor; + m_denom /= factor; + + return *this; + } + + double frac::eval() const + { + double temp(m_num); + return temp /= m_denom; + } + + int frac::num() const + { + return m_num; + } + + unsigned int frac::denom() const + { + return m_denom; + } + + bool frac::is_int() const + { + return __epix_gcd(m_num, m_denom) == m_denom; + } + + + frac operator+ (frac arg1, const frac& arg2) + { + return arg1 += arg2; + } + + frac operator- (frac arg1) + { + return arg1 *= -1; + } + + frac operator- (frac arg1, const frac& arg2) + { + return arg1 -= arg2; + } + + frac operator* (frac arg1, const frac& arg2) + { + return arg1 *= arg2; + } + + frac operator/ (frac arg1, const frac& arg2) + { + return arg1 /= arg2; + } + + // (in)equality + bool operator == (const frac& u, const frac& v) + { + return u.num()*v.denom() == v.num()*u.denom(); + } + + bool operator != (const frac& u, const frac& v) + { + return !(u == v); + } + + // denoms are unsigned + bool operator < (const frac& u, const frac& v) + { + return u.num()*v.denom() < v.num()*u.denom(); + } + + bool operator > (const frac& u, const frac& v) + { + return u.num()*v.denom() > v.num()*u.denom(); + } + + bool operator <= (const frac& u, const frac& v) + { + return u.num()*v.denom() <= v.num()*u.denom(); + } + + bool operator >= (const frac& u, const frac& v) + { + return u.num()*v.denom() >= v.num()*u.denom(); + } + + // N.B.: gcd(0,i) = |i| + unsigned int __epix_gcd (int i, unsigned int j) + { + unsigned int new_i(i<0 ? -i : i); + unsigned int temp; + + if (new_i==0 || j==0) // (1,0) and (0,1) coprime, others not + return new_i + j; + + else { + if (j < new_i) // swap them + { + temp = j; + j = new_i; + new_i = temp; + } + // Euclidean algorithm + while ((temp = j%new_i)) // i does not evenly divide j + { + j = new_i; + new_i = temp; + } + + return new_i; + } + } +} // end of namespace + |