Paranoid Numbers are hard to simplify properly. But they work. Kind of...
# TODO: stb_truetype doesn't compile with some of these warnings.
CXX = g++ -std=gnu++0x -g -Wall -Werror -Wshadow -pedantic -rdynamic
MAIN = main.o
-OBJ = log.o real.o bezier.o document.o objectrenderer.o view.o screen.o vfpu.o quadtree.o graphicsbuffer.o framebuffer.o shaderprogram.o stb_truetype.o gl_core44.o add_digits_asm.o sub_digits_asm.o mul_digits_asm.o div_digits_asm.o arbint.o moc_controlpanel.o controlpanel.o path.o paranoidnumber.o
+OBJ = log.o real.o bezier.o document.o objectrenderer.o view.o screen.o graphicsbuffer.o framebuffer.o shaderprogram.o stb_truetype.o gl_core44.o path.o paranoidnumber.o
QT_INCLUDE := -I/usr/share/qt4/mkspecs/linux-g++-64 -I. -I/usr/include/qt4/QtCore -I/usr/include/qt4/QtGui -I/usr/include/qt4 -I. -Itests -I.
QT_DEF := -DQT_NO_DEBUG -DQT_GUI_LIB -DQT_CORE_LIB
QT_LIB := -L/usr/lib/x86_64-linux-gnu -lQtGui -lQtCore -lpthread
-LIB_x86_64 = ../contrib/lib/libSDL2-2.0.so.0 -lGL -lmpfr -lgmp $(QT_LIB)
-LIB_i386 = ../contrib/lib32/libSDL2-2.0.so.0 -lGL -lgmp
+LIB_x86_64 = ../contrib/lib/libSDL2-2.0.so.0 -lGL
+LIB_i386 = ../contrib/lib32/libSDL2-2.0.so.0 -lGL
LIB_i686 = $(LIB_i386)
MAINRPATH_x86_64 = -Wl,-rpath,'$$ORIGIN/../contrib/lib'
TESTRPATH_i686 = $(TESTRPATH_i386)
OBJPATHS = $(OBJ:%=../obj/%)
DEPS := $(OBJPATHS:%.o=%.d)
-CFLAGS_x86_64 := -I../contrib/include/SDL2 -I`pwd` $(QT_INCLUDE)
+CFLAGS_x86_64 := -I../contrib/include/SDL2 -I`pwd`
CFLAGS_i386 := -I../contrib/include32/SDL2 -I`pwd`
CFLAGS_i686 := $(CFLAGS_i386)
TESTRPATH := $(TESTRPATH_$(ARCH))
CFLAGS := $(CFLAGS_$(ARCH))
-REALTYPE=1
-LINKOBJ = $(OBJPATHS)
-
RM = rm -f
BIN = ../bin/ipdf
+REALTYPE=1
+CONTROLPANEL=enabled
+QUADTREE=disabled
+DEF = -DREALTYPE=$(REALTYPE)
+
+## Only link with things we care about
+
+ifeq ($(QUADTREE),enabled)
+ OBJ := $(OBJ) quadtree.o
+else
+ DEF := $(DEF) -DQUADTREE_DISABLED
+endif
+
+ifeq ($(CONTROLPANEL),enabled)
+ LIB := $(LIB) $(QT_LIB)
+ DEF := $(DEF) $(QT_DEF)
+ CFLAGS := $(CFLAGS) $(QT_INCLUDE)
+ OBJ := $(OBJ) moc_controlpanel.o controlpanel.o
+else
+ DEF := $(DEF) -DCONTROLPANEL_DISABLED
+endif
+
+ifeq ($REALTYPE),3)
+ OBJ := $(OBJ) vfpu.o
+endif
+
+ifeq ($(REALTYPE),5)
+ OBJ := $(OBJ) add_digits_asm.o sub_digits_asm.o mul_digits_asm.o div_digits_asm.o arbint.o
+ LIB := $(LIB) -lgmp
+endif
+
+ifeq ($(REALTYPE),6)
+ LIB := $(LIB) -lgmp -lmpfr
+endif
+
+ifeq ($(REALTYPE),7)
+ LIB := $(LIB) -L../contrib/iRRAM/lib -liRRAM -lgmp -lmpfr
+ CFLAGS := $(CFLAGS) -I../contrib/iRRAM/include
+endif
+
+LINKOBJ = $(OBJPATHS)
-all : REAL = 1
+all : REALTYPE = 1
all : $(BIN)
-single : REAL = 0
+single : REALTYPE = 0
single : $(BIN)
-double : REAL = 1
+double : REALTYPE = 1
double : $(BIN)
-DEF = -DREAL=$(REALTYPE) -DQUADTREE_DISABLED
+
CFLAGS := $(CFLAGS) $(QT_DEF)
$(RM) ../data/movie.ogv
./ipdf | tee ../data/performance.dat | ../tools/stream_plot.py 2>/dev/null & recordmydesktop --fps 10 --on-the-fly-encoding -o ../data/movie.ogv
-# The tests will compile with the default REAL definition
-# To change that you can run as `make DEFS="REAL=X" tests/<target>` where X is your chosen type
+# The tests will compile with the default REALTYPE definition
+# To change that you can run as `make DEFS="REALTYPE=X" tests/<target>` where X is your chosen type
# But remember to make clean first.
tests/% : tests/%.cpp ../obj/tests/%.o $(LINKOBJ)
int main(int argc, char ** argv)
{
+ #if REALTYPE == REAL_IRRAM
+ iRRAM_initialize(argc,argv);
+ #endif
+
#ifndef __STDC_IEC_559__
Warn("__STDC_IEC_559__ not defined. IEEE 754 floating point not fully supported.\n");
#endif
// We want to crash if we ever get a NaN.
feenableexcept(FE_DIVBYZERO | FE_INVALID | FE_OVERFLOW);
- Debug("Compiled with REAL = %d => \"%s\" sizeof(Real) == %d bytes", REAL, g_real_name[REAL], sizeof(Real));
+ Debug("Compiled with REAL = %d => \"%s\" sizeof(Real) == %d bytes", REALTYPE, g_real_name[REALTYPE], sizeof(Real));
Document doc("","fonts/ComicSans.ttf");
srand(time(NULL));
}
if (buttons)
{
- #if REAL == REAL_RATIONAL
+ #if REALTYPE == REAL_RATIONAL
view->Translate(Real(oldx, scr->ViewportWidth()) -Real(x,scr->ViewportWidth()), Real(oldy, scr->ViewportHeight()) - Real(y,scr->ViewportHeight()));
#else
view->Translate(Real(oldx-x)/Real(scr->ViewportWidth()), Real(oldy-y)/Real(scr->ViewportHeight()));
if (wheel)
{
- #if REAL == REAL_RATIONAL
+ #if REALTYPE == REAL_RATIONAL
view->ScaleAroundPoint(Real(x,scr->ViewportWidth()), Real(y,scr->ViewportHeight()), Real(20-wheel, 20));
#else
view->ScaleAroundPoint(Real(x)/Real(scr->ViewportWidth()),Real(y)/Real(scr->ViewportHeight()), Real(expf(-wheel/20.f)));
scr.DebugFontPrintF("[GPU] Render took %lf ms (%lf FPS) (total %lf s, avg FPS %lf)\n", gpu_frame*1e3, 1.0/gpu_frame, total_gpu_time, frames/total_gpu_time);
scr.DebugFontPrintF("[REALTIME] Render+Present+Cruft took %lf ms (%lf FPS) (total %lf s, avg FPS %lf)\n", real_frame*1e3, 1.0/real_frame, total_real_time,frames/total_real_time);
scr.DebugFontPrintF("View bounds: %s\n", view.GetBounds().Str().c_str());
- scr.DebugFontPrintF("type of Real == %s\n", g_real_name[REAL]);
+ scr.DebugFontPrintF("type of Real == %s\n", g_real_name[REALTYPE]);
+ #if REALTYPE == REAL_MPFRCPP
+ scr.DebugFontPrintf("Precision: %s\nRounding: %s\n");
+ #endif
if (view.UsingGPUTransform())
{
{
// hack...
Rect result = view.TransformToViewCoords(Rect(point.x, point.y,1,1));
- int64_t x = Int64(result.x*target.w);
- int64_t y = Int64(result.y*target.h);
+ int64_t x = Int64(result.x)*target.w;
+ int64_t y = Int64(result.y)*target.h;
return PixelPoint(x,y);
}
using namespace std;
namespace IPDF
{
-
+int64_t ParanoidNumber::g_count = 0;
ParanoidNumber::ParanoidNumber(const char * str) : m_value(0), m_op(ADD), m_next_term(NULL), m_next_factor(NULL)
{
//Debug("Now at %s", Str().c_str());
}
- Debug("Constructed {%s} from %s (%f)", Str().c_str(), str, ToDouble());
+ //Debug("Constructed {%s} from %s (%f)", Str().c_str(), str, ToDouble());
}
ParanoidNumber & ParanoidNumber::operator=(const ParanoidNumber & a)
ParanoidNumber & ParanoidNumber::operator+=(const ParanoidNumber & a)
{
- // this = v + t + (a)
- // -> v + (a) + t
+ if (m_next_factor == NULL && a.Floating())
+ {
+ if (ParanoidOp<digit_t>(m_value, a.m_value, ADD))
+ {
+ Simplify();
+ return *this;
+ }
+ }
ParanoidNumber * nt = m_next_term;
ParanoidNumber * nf = m_next_factor;
{
// this = v + t + (a)
// -> v + (a) + t
-
+ if (m_next_factor == NULL && a.Floating())
+ {
+ if (ParanoidOp<digit_t>(m_value, a.m_value, ADD))
+ {
+ Simplify();
+ return *this;
+ }
+ }
+
ParanoidNumber * nt = m_next_term;
ParanoidNumber * nf = m_next_factor;
ParanoidNumber & ParanoidNumber::operator*=(const ParanoidNumber & a)
{
- Debug("{%s} *= {%s}", Str().c_str(), a.Str().c_str());
+
+ //if (m_value == 0)
+ // return *this;
+ //Debug("{%s} *= {%s}", Str().c_str(), a.Str().c_str());
// this = (vf + t) * (a)
- ParanoidNumber * nf = m_next_factor;
- m_next_factor = new ParanoidNumber(a, MULTIPLY);
- ParanoidNumber * t = m_next_factor;
+ if (a.Floating() && ParanoidOp<digit_t>(m_value, a.m_value, MULTIPLY))
+ {
+ if (m_next_term != NULL)
+ m_next_term->operator*=(a);
+ Simplify();
+ return *this;
+ }
+
+ ParanoidNumber * t = this;
while (t->m_next_factor != NULL)
t = t->m_next_factor;
- t->m_next_factor = nf;
+ t->m_next_factor = new ParanoidNumber(a, MULTIPLY);
+
if (m_next_term != NULL)
m_next_term->operator*=(a);
- //Debug("Simplify after mul");
- Debug("Simplify {%s}", Str().c_str());
+
+ //Debug("Simplify {%s}", Str().c_str());
Simplify();
+ //Debug("Simplified to {%s}", Str().c_str());
return *this;
}
ParanoidNumber & ParanoidNumber::operator/=(const ParanoidNumber & a)
{
+
+
+
+ if (a.Floating() && ParanoidOp<digit_t>(m_value, a.m_value, DIVIDE))
+ {
+ if (m_next_term != NULL)
+ m_next_term->operator/=(a);
+ Simplify();
+ return *this;
+ }
+
//Debug("Called %s /= %s", Str().c_str(), a.Str().c_str());
// this = (vf + t) * (a)
- ParanoidNumber * nf = m_next_factor;
- m_next_factor = new ParanoidNumber(a, DIVIDE);
- ParanoidNumber * t = m_next_factor;
+ ParanoidNumber * t = this;
while (t->m_next_factor != NULL)
+ {
t = t->m_next_factor;
- t->m_next_factor = nf;
+ }
+ t->m_next_factor = new ParanoidNumber(a, DIVIDE);
+
if (m_next_term != NULL)
m_next_term->operator/=(a);
+
Simplify();
return *this;
}
void ParanoidNumber::SimplifyTerms()
{
+
//Debug("Simplify {%s}", Str().c_str());
if (m_next_term == NULL)
{
for (ParanoidNumber * a = this; a != NULL; a = a->m_next_term)
{
- ParanoidNumber * bprev = a;
ParanoidNumber * b = a->m_next_term;
+ if (a->m_next_factor != NULL)
+ {
+ continue;
+ }
+
+ ParanoidNumber * bprev = a;
while (b != NULL)
{
//Debug("Simplify factors of %s", b->Str().c_str());
b = b->m_next_term;
continue;
}
- float f = a->m_value;
- feclearexcept(FE_ALL_EXCEPT);
- switch (b->m_op)
- {
- case ADD:
- f += b->m_value;
- break;
- case SUBTRACT:
- f -= b->m_value;
- break;
- default:
- Fatal("Unexpected %c in term list...", OpChar(b->m_op));
- break;
- }
- if (!fetestexcept(FE_ALL_EXCEPT))
+ bool simplify = false;
+ simplify = ParanoidOp<digit_t>(a->m_value, b->Head<digit_t>(), ADD);
+ if (simplify)
{
- a->m_value = f;
bprev->m_next_term = b->m_next_term;
b->m_next_term = NULL;
delete b;
b = bprev;
}
+
bprev = b;
b = b->m_next_term;
}
void ParanoidNumber::SimplifyFactors()
{
+
//Debug("Simplify {%s}", Str().c_str());
if (m_next_factor == NULL)
{
for (ParanoidNumber * a = this; a != NULL; a = a->m_next_factor)
{
+ if ((a->m_op != ADD || a->m_op != SUBTRACT) && a->m_next_term != NULL)
+ continue;
+
ParanoidNumber * bprev = a;
ParanoidNumber * b = a->m_next_factor;
while (b != NULL)
b = b->m_next_factor;
continue;
}
- float f = a->m_value;
- feclearexcept(FE_ALL_EXCEPT);
- switch (b->m_op)
+
+ Optype op = b->m_op;
+ if (a->m_op == DIVIDE)
{
- case MULTIPLY:
- if (a->m_op != DIVIDE)
- f *= b->m_value;
- else
- f /= b->m_value;
- break;
- case DIVIDE:
- if (a->m_op != DIVIDE)
- f /= b->m_value;
- else
- f *= b->m_value;
- break;
- default:
- Fatal("Unexpected %c in factor list...",OpChar(b->m_op));
- break;
+ op = (b->m_op == DIVIDE) ? MULTIPLY : DIVIDE;
}
- if (!fetestexcept(FE_ALL_EXCEPT))
- {
-
- a->m_value = f;
+
+ if (ParanoidOp<digit_t>(a->m_value, b->m_value, op))
+ {
+
bprev->m_next_factor = b->m_next_factor;
b->m_next_factor = NULL;
delete b;
b = bprev;
}
- //else
- //Debug("Failed to simplify %f %c %f", a->m_value, OpChar(b->m_op), b->m_value);
bprev = b;
b = b->m_next_factor;
}
{
string result("");
stringstream s;
-
- s << m_value;
+ s << (double)m_value;
if (m_next_factor != NULL)
{
- result += OpChar(m_op);
- result += "(";
result += s.str();
- result += m_next_factor->Str();
- result += ")";
+ result += OpChar(m_next_factor->m_op);
+ if (m_next_factor->m_next_term != NULL)
+ result += "(" + m_next_factor->Str() + ")";
+ else
+ result += m_next_factor->Str();
}
else
{
- result += OpChar(m_op);
result += s.str();
}
if (m_next_term != NULL)
{
result += " ";
+ result += OpChar(m_next_term->m_op);
result += m_next_term->Str();
}
return result;
}
+template <>
+bool TrustingOp<float>(float & a, const float & b, Optype op)
+{
+ feclearexcept(FE_ALL_EXCEPT);
+ switch (op)
+ {
+ case ADD:
+ a += b;
+ break;
+ case SUBTRACT:
+ a -= b;
+ break;
+ case MULTIPLY:
+ a *= b;
+ break;
+ case DIVIDE:
+ a /= b;
+ break;
+ }
+ return !fetestexcept(FE_ALL_EXCEPT);
+}
+
+template <>
+bool TrustingOp<double>(double & a, const double & b, Optype op)
+{
+ feclearexcept(FE_ALL_EXCEPT);
+ switch (op)
+ {
+ case ADD:
+ a += b;
+ break;
+ case SUBTRACT:
+ a -= b;
+ break;
+ case MULTIPLY:
+ a *= b;
+ break;
+ case DIVIDE:
+ a /= b;
+ break;
+ }
+ return !fetestexcept(FE_ALL_EXCEPT);
+}
+
+template <>
+bool TrustingOp<int8_t>(int8_t & a, const int8_t & b, Optype op)
+{
+ int16_t sa(a);
+ bool exact = true;
+ switch (op)
+ {
+ case ADD:
+ sa += b;
+ exact = (abs(sa) <= 127);
+ break;
+ case SUBTRACT:
+ sa -= b;
+ exact = (abs(sa) <= 127);
+ break;
+ case MULTIPLY:
+ sa *= b;
+ exact = (abs(sa) <= 127);
+ break;
+ case DIVIDE:
+ exact = (b != 0 && sa > b && sa % b == 0);
+ sa /= b;
+ break;
+ }
+ a = (int8_t)(sa);
+ return exact;
+}
+
}
#include <map>
#include <string>
#include "log.h"
+#include <fenv.h>
-
+#define PARANOID_DIGIT_T int8_t // we could theoretically replace this with a template
+ // but let's not do that...
namespace IPDF
{
+ typedef enum {ADD, SUBTRACT, MULTIPLY, DIVIDE} Optype;
+
+ /** Performs an operation, returning if the result was exact **/
+ // NOTE: DIFFERENT to ParanoidOp (although that wraps to this...)
+ template <class T> bool TrustingOp(T & a, const T & b, Optype op);
+
+ /** Performs an operation _only_ if the result would be exact **/
+ template <class T> bool ParanoidOp(T & a, const T & b, Optype op)
+ {
+ T cpy(a);
+ if (TrustingOp<T>(cpy, b, op))
+ {
+ a = cpy;
+ return true;
+ }
+ return false;
+ }
+
+
+ template <> bool TrustingOp<float>(float & a, const float & b, Optype op);
+ template <> bool TrustingOp<double>(double & a, const double & b, Optype op);
+ template <> bool TrustingOp<int8_t>(int8_t & a, const int8_t & b, Optype op);
+
+ // Attempt to comine two terms: a*b + c*d or a/b + c/d
+ template <class T> bool CombineTerms(T & aa, Optype aop, T & bb, T & cc, Optype cop, T & dd)
+ {
+ T a(aa); T b(bb); T c(cc); T d(dd);
+ if (aop == MULTIPLY && cop == MULTIPLY) // a*b + c*d
+ {
+ if ((ParanoidOp<T>(c, b, DIVIDE) || ParanoidOp(d, b, DIVIDE))
+ && TrustingOp<T>(c, d, MULTIPLY) && TrustingOp<T>(a,c,ADD)
+ && TrustingOp<T>(a, b, MULTIPLY)) // (a + (cd)/b) * b
+ {
+ aa = a;
+ bb = 1;
+ cc = 1;
+ dd = 1;
+ return true;
+ }
+ if ((ParanoidOp<T>(a, d, DIVIDE) || ParanoidOp(b, d, DIVIDE))
+ && TrustingOp<T>(a, b, MULTIPLY) && TrustingOp<T>(a,c,ADD)
+ && TrustingOp<T>(a, d, MULTIPLY)) // ((ab)/d + c)*d
+ {
+ aa = a;
+ bb = 1;
+ cc = 1;
+ dd = 1;
+ return true;
+ }
+ return false;
+ }
+ else if (aop == DIVIDE && cop == DIVIDE)
+ {
+ if (TrustingOp<T>(a, d, MULTIPLY) && TrustingOp<T>(c, b, MULTIPLY)
+ && TrustingOp<T>(a, c, ADD) && TrustingOp<T>(b, d, MULTIPLY))
+ {
+ cc = 1;
+ dd = 1;
+ if (ParanoidOp<T>(a, b, DIVIDE))
+ {
+ aa = a;
+ bb = 1;
+ return true;
+ }
+ aa = a;
+ bb = b;
+ return true;
+ }
+ return false;
+ }
+ return false;
+ }
+
class ParanoidNumber
{
+
public:
- typedef enum {ADD, SUBTRACT, MULTIPLY, DIVIDE} Optype;
-
- ParanoidNumber(float value=0, Optype type = ADD) : m_value(value), m_op(type), m_next_term(NULL), m_next_factor(NULL)
+ typedef PARANOID_DIGIT_T digit_t;
+
+ ParanoidNumber(digit_t value=0, Optype type = ADD) : m_value(value), m_op(type), m_next_term(NULL), m_next_factor(NULL)
{
-
+ Construct();
}
ParanoidNumber(const ParanoidNumber & cpy) : m_value(cpy.m_value), m_op(cpy.m_op), m_next_term(NULL), m_next_factor(NULL)
{
m_next_factor = new ParanoidNumber(*(cpy.m_next_factor));
}
+ Construct();
}
- ParanoidNumber(const ParanoidNumber & cpy, Optype type) : m_value(cpy.m_value), m_op(type), m_next_term(NULL), m_next_factor(NULL)
+ ParanoidNumber(const ParanoidNumber & cpy, Optype type) : ParanoidNumber(cpy)
{
- if (cpy.m_next_term != NULL)
- {
- m_next_term = new ParanoidNumber(*(cpy.m_next_term));
- }
- if (cpy.m_next_factor != NULL)
- {
- m_next_factor = new ParanoidNumber(*(cpy.m_next_factor));
- }
+ m_op = type;
}
ParanoidNumber(const char * str);
- ParanoidNumber(const std::string & str) : ParanoidNumber(str.c_str()) {}
+ ParanoidNumber(const std::string & str) : ParanoidNumber(str.c_str()) {Construct();}
virtual ~ParanoidNumber()
{
delete m_next_term;
if (m_next_factor != NULL)
delete m_next_factor;
+ g_count--;
}
+ inline void Construct() {g_count++;}
+
+
template <class T> T Convert() const;
+ template <class T> T AddTerms() const;
+ template <class T> T MultiplyFactors() const;
+ template <class T> T Head() const {return (m_op == SUBTRACT) ? T(-m_value) : T(m_value);}
+
+
+
+
double ToDouble() const {return Convert<double>();}
float ToFloat() const {return Convert<float>();}
+ digit_t Digit() const {return Convert<digit_t>();}
bool Floating() const {return (m_next_term == NULL && m_next_factor == NULL);}
bool Sunken() const {return !Floating();} // I could not resist...
return opch[(int)op];
}
+ static int64_t Paranoia() {return g_count;}
+
private:
+ static int64_t g_count;
void Simplify();
void SimplifyTerms();
void SimplifyFactors();
- float m_value;
+ digit_t m_value;
Optype m_op;
ParanoidNumber * m_next_term;
ParanoidNumber * m_next_factor;
-
-
-
-
-
-
};
-
template <class T>
-T ParanoidNumber::Convert() const
+T ParanoidNumber::AddTerms() const
{
- T value = (m_op == SUBTRACT) ? -m_value : m_value;
- const ParanoidNumber * n = m_next_factor;
- if (n != NULL)
+ T value(0);
+ for (ParanoidNumber * a = m_next_term; a != NULL; a = a->m_next_term)
{
- switch (n->m_op)
- {
- case MULTIPLY:
- value *= n->Convert<T>();
- break;
- case DIVIDE:
- value /= n->Convert<T>();
- break;
- default:
- Fatal("Shouldn't happen");
- break;
- }
+ value += a->Head<T>() * a->MultiplyFactors<T>();
}
- n = m_next_term;
- if (n != NULL)
+ return value;
+}
+
+template <class T>
+T ParanoidNumber::MultiplyFactors() const
+{
+ T value(1);
+ for (ParanoidNumber * a = m_next_factor; a != NULL; a = a->m_next_factor)
{
- switch (n->m_op)
- {
- case ADD:
- case SUBTRACT:
- value += n->Convert<T>();
- break;
- default:
- Fatal("Shouldn't happen");
- }
+ if (a->m_op == DIVIDE)
+ value /= (a->Head<T>() + a->AddTerms<T>());
+ else
+ value *= (a->Head<T>() + a->AddTerms<T>());
}
return value;
}
+
+template <class T>
+T ParanoidNumber::Convert() const
+{
+ return Head<T>() * MultiplyFactors<T>() + AddTerms<T>();
+}
+
+
+
}
#endif //_PARANOIDNUMBER_H
--- /dev/null
+#include <sys/time.h>
+#include <unistd.h>
+#include <cstdio>
+
+/**
+ * Print a progress bar to stderr + elapsed time since execution
+ * Stolen from some of my Codejam stuff, but modified to actually print the correct times...
+ */
+inline void ProgressBar(int64_t c, int64_t nCases, int len=20)
+{
+ static timeval startTime;
+ static timeval thisTime;
+ static timeval lastTime;
+ static int64_t total_usec = 0;
+ gettimeofday(&thisTime, NULL);
+ if (c <= 0)
+ {
+ startTime = thisTime;
+ lastTime = thisTime;
+ total_usec = 0;
+ }
+
+
+ fprintf(stderr, "\33[2K\r");
+ //fprintf(stderr, "\n");
+
+ fprintf(stderr, "[");
+ int i = 0;
+ for (i = 0; i < ((len*c)/nCases); ++i)
+ fprintf(stderr, "=");
+ fprintf(stderr, ">");
+ for (; i < len; ++i)
+ fprintf(stderr, " ");
+
+
+ int64_t delta_usec = 1000000*(thisTime.tv_sec-lastTime.tv_sec) + (thisTime.tv_usec-lastTime.tv_usec);
+ total_usec += delta_usec;
+
+ int64_t ds = delta_usec / 1000000;
+ int64_t dus = delta_usec - 1000000*ds;
+ int64_t ts = total_usec / 1000000;
+ int64_t tus = total_usec - 1000000*ts;
+
+ fprintf(stderr, "] Case #%li: %02li:%02li.%06li Total: %02li:%02li.%06li| ", c,
+ ds/60, ds%60, dus,
+ ts/60, ts%60, tus);
+
+ lastTime = thisTime;
+}
#include "real.h"
+#include <fenv.h>
namespace IPDF
{
- // Maps the REAL to a string
+ // Maps the REALTYPE to a string
const char * g_real_name[] = {
"single",
"double",
"VFPU",
"Rational<int64_t>",
"Rational<Arbint>",
- "mpfrc++ real"
+ "mpfrc++ real",
+ "iRRAM REAL"
};
-#if REAL == REAL_RATIONAL_ARBINT
+#if REALTYPE == REAL_RATIONAL_ARBINT
template <> Gmpint Tabs(const Gmpint & a)
{
return a.Abs();
}
#endif
+
+
+
}
#define REAL_RATIONAL 4
#define REAL_RATIONAL_ARBINT 5
#define REAL_MPFRCPP 6
+#define REAL_IRRAM 7
-#ifndef REAL
- #error "REAL was not defined!"
+#ifndef REALTYPE
+ #error "REALTYPE was not defined!"
#endif
-#if REAL == REAL_VFPU
+#if REALTYPE == REAL_VFPU
#include "vfpu.h"
#endif
-#if REAL == REAL_RATIONAL
+#if REALTYPE == REAL_RATIONAL
#include "rational.h"
-#endif //REAL
+#endif //REALTYPE
-#if REAL == REAL_RATIONAL_ARBINT
+#if REALTYPE == REAL_RATIONAL_ARBINT
#include "rational.h"
#include "arbint.h"
#include "gmpint.h"
-#endif //REAL
+#endif //REALTYPE
-#if REAL == REAL_MPFRCPP
+#if REALTYPE == REAL_MPFRCPP
#include <mpreal.h>
-#endif //REAL
+#endif //REALTYPE
+
+#if REALTYPE == REAL_IRRAM
+ #include "../contrib/iRRAM/include/iRRAM/lib.h"
+#endif
namespace IPDF
{
extern const char * g_real_name[];
-#if REAL == REAL_SINGLE
+#if REALTYPE == REAL_SINGLE
typedef float Real;
-#elif REAL == REAL_DOUBLE
+#elif REALTYPE == REAL_DOUBLE
typedef double Real;
-#elif REAL == REAL_LONG_DOUBLE
+#elif REALTYPE == REAL_LONG_DOUBLE
typedef long double Real;
-#elif REAL == REAL_VFPU
+#elif REALTYPE == REAL_VFPU
typedef VFPU::VFloat Real;
inline float Float(const Real & r) {return r.m_value;}
inline double Double(const Real & r) {return r.m_value;}
-#elif REAL == REAL_RATIONAL
+#elif REALTYPE == REAL_RATIONAL
typedef Rational<int64_t> Real;
inline float Float(const Real & r) {return (float)r.ToDouble();}
inline double Double(const Real & r) {return r.ToDouble();}
-#elif REAL == REAL_RATIONAL_ARBINT
+#elif REALTYPE == REAL_RATIONAL_ARBINT
#define ARBINT Gmpint // Set to Gmpint or Arbint here
typedef Rational<ARBINT> Real;
inline double Double(const Real & r) {return r.ToDouble();}
inline int64_t Int64(const Real & r) {return r.ToInt64();}
inline Rational<ARBINT> Sqrt(const Rational<ARBINT> & r) {return r.Sqrt();}
-#elif REAL == REAL_MPFRCPP
+#elif REALTYPE == REAL_MPFRCPP
typedef mpfr::mpreal Real;
inline double Double(const Real & r) {return r.toDouble();}
inline float Float(const Real & r) {return r.toDouble();}
inline int64_t Int64(const Real & r) {return r.toLong();}
inline Real Sqrt(const Real & r) {return mpfr::sqrt(r, mpfr::mpreal::get_default_rnd());}
inline Real Abs(const Real & r) {return mpfr::abs(r, mpfr::mpreal::get_default_rnd());}
+#elif REALTYPE == REAL_IRRAM
+ typedef iRRAM::REAL Real;
+ inline double Double(const Real & r) {return r.as_double(53);}
+ inline float Float(const Real & r) {return r.as_double(53);}
+ inline int64_t Int64(const Real & r) {return (int64_t)r.as_double(53);}
+ inline Real Sqrt(const Real & r) {return iRRAM::sqrt(r);}
#else
#error "Type of Real unspecified."
-#endif //REAL
+#endif //REALTYPE
// Allow us to call Float on the primative types
// Useful so I can template some things that could be either (a more complicated) Real or a primitive type
// Mostly in the testers.
- inline float Float(float f) {return (float)f;}
inline float Float(double f) {return (float)f;}
inline float Float(long double f) {return (float)(f);}
inline double Double(float f) {return (double)f;}
Vec2() : x(0), y(0) {}
Vec2(Real _x, Real _y) : x(_x), y(_y) {}
Vec2(const std::pair<Real, Real> & p) : x(p.first), y(p.second) {}
+ #if REALTYPE != REAL_IRRAM
Vec2(const std::pair<int64_t, int64_t> & p) : x(p.first), y(p.second) {}
-
+ #else
+ Vec2(const std::pair<int64_t, int64_t> & p) : x((int)p.first), y((int)p.second) {}
+ // Apparently iRRAM didn't implement -= for a constant argument
+ #endif
bool operator==(const Vec2& other) const { return (x == other.x) && (y == other.y); }
bool operator!=(const Vec2& other) const { return !(*this == other); }
};
+
+
}
#endif //_REAL_H
int main(int argc, char ** argv)
{
+ #if REALTYPE == REAL_IRRAM
+ iRRAM_initialize(argc,argv);
+ #endif
while (cin.good())
{
double da; double db;
Debug("FLT_MAX = %.40f", FLT_MAX);
Debug("FLT_MIN = %.40f", FLT_MIN);
Debug("FLT_EPSILON = %.40f", FLT_EPSILON);
-
- ParanoidNumber a("0.3");
- Debug("start at {%s} = %lf", a.Str().c_str(), a.ToDouble());
- cout << "0.3 ";
- float fa = 0.3;
- double da = 0.3;
+ Debug("Sizeof ParanoidNumber::digit_t is %u", sizeof(ParanoidNumber::digit_t));
+ Debug("Sizeof ParanoidNumber is %u", sizeof(ParanoidNumber));
+
+ string token("");
+ cin >> token;
+ ParanoidNumber a(token.c_str());
+ double da = a.ToDouble();
+ float fa = da;
while (cin.good())
{
char op;
cin >> op;
- string token("");
+ token = "";
for (char c = cin.peek(); cin.good() && !iswspace(c); c = cin.peek())
{
if (c == '+' || c == '-' || c == '*' || c == '/')
#include <cfloat>
#include <fenv.h>
#include "paranoidnumber.h"
+#include "progressbar.h"
using namespace std;
using namespace IPDF;
-string RandomNumberAsString(int digits = 6)
+string RandomNumberAsString(int max_digits = 12)
{
string result("");
- int dp = 1+(rand() % 3);
+ int digits = 1+(rand() % max_digits);
+ int dp = (rand() % digits)+1;
for (int i = 0; i < digits; ++i)
{
if (i == dp)
return result;
}
-#define TEST_CASES 10000
+bool CloseEnough(double d, ParanoidNumber & p)
+{
+ double pd = p.ToDouble();
+
+ if (d == 0)
+ return fabs(pd) <= 1e-6;
+ return fabs((fabs(pd - d) / d)) <= 1e-6;
+}
+
+#define TEST_CASES 1000
int main(int argc, char ** argv)
{
-
+ srand(time(NULL));
string number(RandomNumberAsString());
ParanoidNumber a(number);
float fa = strtof(number.c_str(), NULL);
double da = strtod(number.c_str(), NULL);
+ double diff = 0;
long double lda = strtold(number.c_str(), NULL);
if (fabs(a.ToDouble() - da) > 1e-6)
}
char opch[] = {'+','-','*','/'};
-
+ int opcount[] = {0,0,0,0};
for (int i = 0; i < TEST_CASES; ++i)
{
+ ProgressBar(i, TEST_CASES); fprintf(stderr, "%.30lf (%d)+ (%d)- (%d)* (%d)/ [Paranoia %ld]",diff, opcount[0],opcount[1],opcount[2],opcount[3], ParanoidNumber::Paranoia());
number = RandomNumberAsString();
ParanoidNumber b(number);
float fb = strtof(number.c_str(), NULL);
double db = strtod(number.c_str(), NULL);
+ if (db == 0)
+ {
+ --i;
+ continue;
+ }
long double ldb = strtold(number.c_str(), NULL);
- int op = (rand() % 4);
+ int op = rand() % 4;//= 2*(rand() % 2);
+ while (op == 1)
+ op = rand() % 4;
+ //if (rand() % 2 == 0)
+ //op = 2+(rand() % 2);
+
+ opcount[op]++;
ParanoidNumber olda(a);
double oldda(da);
switch (op)
lda /= ldb;
break;
}
- if (fabs(a.ToDouble() - da) > 1.0 )
+ diff = 100.0*(fabs(a.ToDouble() - da) / da);
+ if (!CloseEnough(da, a))
{
Error("Op %i: ParanoidNumber probably doesn't work", i);
Error("Operation: %lf %c %lf", oldda, opch[op], db);
Error("As PN: %lf %c %lf", olda.ToDouble(), opch[op], b.ToDouble());
- Fatal("%lf, expected aboout %lf", a.ToDouble(), da);
+ Error("PN String: %s", a.Str().c_str());
+ Error("Diff is %.40lf", diff);
+ Fatal("%.40lf, expected aboout %.40lf", a.ToDouble(), da);
+
+
}
+
}
printf("ParanoidNumber: {%s} = %.40lf\n", a.Str().c_str(), a.ToDouble());
printf("float: %.40f\n", fa);
printf("double: %.40lf\n", da);
printf("long double: %.40Lf\n", lda);
+ printf("diff %.40lf\n", diff);
}
git.ucc.asn.au / ipdf / code.git / commitdiff
UCC git Repository :: git.ucc.asn.au