X-Git-Url: https://git.ucc.asn.au/?p=ipdf%2Fcode.git;a=blobdiff_plain;f=src%2Fparanoidnumber.cpp;h=60802eb432a7acc0f17c530a42bf7b0461c17f0b;hp=b45c71e069e7ca95b0d79d7b1c8987fc9403080f;hb=64b7c42c71c35d520424cf4ca5ecaa99faef8b26;hpb=77137590512d969da2d54d9ba53d76836a290c6a diff --git a/src/paranoidnumber.cpp b/src/paranoidnumber.cpp index b45c71e..60802eb 100644 --- a/src/paranoidnumber.cpp +++ b/src/paranoidnumber.cpp @@ -6,15 +6,19 @@ #include #include +// here be many copy paste bugs + using namespace std; namespace IPDF { -int64_t ParanoidNumber::g_count = 0; +#ifdef PARANOID_USE_ARENA +ParanoidNumber::Arena ParanoidNumber::g_arena; +#endif //PARANOID_USE_ARENA + ParanoidNumber::~ParanoidNumber() { - g_count--; for (int i = 0; i < NOP; ++i) { for (auto n : m_next[i]) @@ -22,11 +26,23 @@ ParanoidNumber::~ParanoidNumber() } } -ParanoidNumber::ParanoidNumber(const char * str) : m_value(0), m_cached_result(0) +ParanoidNumber::ParanoidNumber(const string & str) : m_value(0), m_next() { - Construct(); + #ifdef PARANOID_SIZE_LIMIT + m_size = 1; + #endif + #ifdef PARANOID_CACHE_RESULTS + m_cache_valid = false; + #endif + int dp = 0; int end = 0; + bool negate = str[0] == '-'; + if (negate) + { + dp++; + end++; + } while (str[dp] != '\0' && str[dp] != '.') { ++dp; @@ -35,7 +51,7 @@ ParanoidNumber::ParanoidNumber(const char * str) : m_value(0), m_cached_result(0 while (str[end] != '\0') ++end; ParanoidNumber m(1); - for (int i = dp-1; i >= 0; --i) + for (int i = dp-1; i >= negate; --i) { ParanoidNumber b(str[i]-'0'); b*=m; @@ -50,31 +66,48 @@ ParanoidNumber::ParanoidNumber(const char * str) : m_value(0), m_cached_result(0 b*=n; this->operator+=(b); } + + if (negate) + Negate(); + + #ifdef PARANOID_COMPARE_EPSILON + double d = strtod(str.c_str(), NULL); + CompareForSanity(d, d); + #endif } ParanoidNumber & ParanoidNumber::operator=(const ParanoidNumber & a) { + //assert(this != NULL); + + #ifdef PARANOID_SIZE_LIMIT + m_size = a.m_size; + #endif + m_value = a.m_value; + #ifdef PARANOID_CACHE_RESULTS m_cached_result = a.m_cached_result; + m_cache_valid = a.m_cache_valid; + #endif for (int i = 0; i < NOP; ++i) { - for (unsigned j = 0; j < m_next[i].size() && j < a.m_next[i].size(); ++j) - { - m_next[i][j]->operator=(*(a.m_next[i][j])); - } - - for (unsigned j = a.m_next[i].size(); j < m_next[i].size(); ++j) - { - delete m_next[i][j]; - } - m_next[i].resize(a.m_next[i].size()); - } + for (auto n : m_next[i]) + delete n; + m_next[i].clear(); + for (auto n : a.m_next[i]) + m_next[i].push_back(new ParanoidNumber(*n)); + } + #ifdef PARANOID_COMPARE_EPSILON + CompareForSanity(a.Digit(),a.Digit()); + #endif return *this; } string ParanoidNumber::Str() const { + + //assert(this != NULL); string result(""); stringstream s; s << (double)m_value; @@ -120,6 +153,8 @@ string ParanoidNumber::Str() const template <> bool TrustingOp(float & a, const float & b, Optype op) { + + feclearexcept(FE_ALL_EXCEPT); switch (op) { @@ -133,6 +168,12 @@ bool TrustingOp(float & a, const float & b, Optype op) a *= b; break; case DIVIDE: + if (b == 0) + { + a = (a >= 0) ? INFINITY : -INFINITY; + return false; + } + a /= b; break; case NOP: @@ -157,6 +198,45 @@ bool TrustingOp(double & a, const double & b, Optype op) a *= b; break; case DIVIDE: + if (b == 0) + { + a = (a >= 0) ? INFINITY : -INFINITY; + return false; + } + a /= b; + break; + case NOP: + break; + } + return !fetestexcept(FE_ALL_EXCEPT); +} + + + +template <> +bool TrustingOp(long double & a, const long 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: + if (b == 0) + { + a = (a >= 0) ? INFINITY : -INFINITY; + return false; + } + a /= b; break; case NOP: @@ -165,6 +245,7 @@ bool TrustingOp(double & a, const double & b, Optype op) return !fetestexcept(FE_ALL_EXCEPT); } + template <> bool TrustingOp(int8_t & a, const int8_t & b, Optype op) { @@ -196,40 +277,179 @@ bool TrustingOp(int8_t & a, const int8_t & b, Optype op) } -ParanoidNumber & ParanoidNumber::operator+=(const ParanoidNumber & a) +ParanoidNumber & ParanoidNumber::operator+=(const digit_t & a) { + #ifdef PARANOID_COMPARE_EPSILON + digit_t compare = Digit(); + compare += a; + #endif delete Operation(new ParanoidNumber(a), ADD); + Simplify(SUBTRACT); + Simplify(ADD); + #ifdef PARANOID_COMPARE_EPSILON + CompareForSanity(compare, a); + #endif + return *this; +} + + +ParanoidNumber & ParanoidNumber::operator-=(const digit_t & a) +{ + #ifdef PARANOID_COMPARE_EPSILON + digit_t compare = Digit(); + compare -= a; + #endif + delete Operation(new ParanoidNumber(a), SUBTRACT); Simplify(ADD); Simplify(SUBTRACT); + #ifdef PARANOID_COMPARE_EPSILON + CompareForSanity(compare, a); + #endif + return *this; +} + +ParanoidNumber & ParanoidNumber::operator*=(const digit_t & a) +{ + #ifdef PARANOID_COMPARE_EPSILON + digit_t compare = Digit(); + compare *= a; + #endif + delete Operation(new ParanoidNumber(a), MULTIPLY); + Simplify(DIVIDE); + Simplify(MULTIPLY); + #ifdef PARANOID_COMPARE_EPSILON + CompareForSanity(compare, a); + #endif + return *this; +} + + +ParanoidNumber & ParanoidNumber::operator/=(const digit_t & a) +{ + #ifdef PARANOID_COMPARE_EPSILON + digit_t compare = Digit(); + compare /= a; + #endif + delete Operation(new ParanoidNumber(a), DIVIDE); + Simplify(MULTIPLY); + Simplify(DIVIDE); + #ifdef PARANOID_COMPARE_EPSILON + CompareForSanity(compare, a); + #endif + return *this; +} + + +ParanoidNumber & ParanoidNumber::operator+=(const ParanoidNumber & a) +{ + #ifdef PARANOID_COMPARE_EPSILON + digit_t compare = Digit(); + compare += a.Digit(); + #endif + delete Operation(new ParanoidNumber(a), ADD); + Simplify(SUBTRACT); + Simplify(ADD); + #ifdef PARANOID_COMPARE_EPSILON + CompareForSanity(compare, a.Digit()); + #endif return *this; } ParanoidNumber & ParanoidNumber::operator-=(const ParanoidNumber & a) { + #ifdef PARANOID_COMPARE_EPSILON + digit_t compare = Digit(); + compare -= a.Digit(); + #endif delete Operation(new ParanoidNumber(a), SUBTRACT); - //Simplify(SUBTRACT); - //Simplify(ADD); + Simplify(ADD); + Simplify(SUBTRACT); + #ifdef PARANOID_COMPARE_EPSILON + CompareForSanity(compare, a.Digit()); + #endif return *this; } ParanoidNumber & ParanoidNumber::operator*=(const ParanoidNumber & a) { + #ifdef PARANOID_COMPARE_EPSILON + digit_t compare = Digit(); + compare *= a.Digit(); + #endif delete Operation(new ParanoidNumber(a), MULTIPLY); + Simplify(DIVIDE); + Simplify(MULTIPLY); + #ifdef PARANOID_COMPARE_EPSILON + CompareForSanity(compare, a.Digit()); + #endif return *this; } ParanoidNumber & ParanoidNumber::operator/=(const ParanoidNumber & a) { + #ifdef PARANOID_COMPARE_EPSILON + digit_t compare = Digit(); + compare /= a.Digit(); + #endif delete Operation(new ParanoidNumber(a), DIVIDE); + Simplify(MULTIPLY); + Simplify(DIVIDE); + #ifdef PARANOID_COMPARE_EPSILON + CompareForSanity(compare, a.Digit()); + #endif + return *this; +} + +ParanoidNumber & ParanoidNumber::operator=(const digit_t & a) +{ + + for (int i = 0; i < NOP; ++i) + { + for (auto n : m_next[i]) + delete n; + m_next[i].clear(); + } + m_value = a; + #ifdef PARANOID_CACHE_RESULTS + m_cached_result = a; + m_cache_valid = true; + #endif + + #ifdef PARANOID_COMPARE_EPSILON + CompareForSanity(a,a); + #endif + return *this; } // a + b ParanoidNumber * ParanoidNumber::OperationTerm(ParanoidNumber * b, Optype op, ParanoidNumber ** merge_point, Optype * merge_op) { - m_cached_result = nan(""); + ////assert(b->SanityCheck()); + #ifdef PARANOID_CACHE_RESULTS + m_cache_valid = false; + #endif + #ifdef PARANOID_SIZE_LIMIT + if (m_size + b->m_size >= PARANOID_SIZE_LIMIT) + { + this->operator=(this->Digit()); + if (op == ADD) + m_value += b->Digit(); + else + m_value -= b->Digit(); + m_size = 1; + #ifdef PARANOID_CACHE_RESULTS + m_cached_result = m_value; + m_cache_valid = true; + #endif + return b; + } + //Debug("At size limit %d", m_size); + #endif + + if (Floating() && m_value == 0) // 0 + b = b { m_value = b->m_value; @@ -244,6 +464,8 @@ ParanoidNumber * ParanoidNumber::OperationTerm(ParanoidNumber * b, Optype op, Pa m_next[i] = b->m_next[i]; b->m_next[i].clear(); } + + //assert(SanityCheck()); return b; } if (b->Floating() && b->m_value == 0) // a + 0 = a @@ -259,14 +481,15 @@ ParanoidNumber * ParanoidNumber::OperationTerm(ParanoidNumber * b, Optype op, Pa Optype addop = (op == ADD) ? ADD : SUBTRACT; for (auto add : b->m_next[ADD]) { - delete OperationTerm(add, addop); + delete (OperationTerm(add, addop)); } Optype subop = (op == ADD) ? SUBTRACT : ADD; for (auto sub : b->m_next[SUBTRACT]) - delete OperationTerm(sub, subop); + delete (OperationTerm(sub, subop)); b->m_next[ADD].clear(); b->m_next[SUBTRACT].clear(); + //assert(SanityCheck()); return b; } } @@ -276,7 +499,7 @@ ParanoidNumber * ParanoidNumber::OperationTerm(ParanoidNumber * b, Optype op, Pa bool parent = (merge_point == NULL); ParanoidNumber * merge = this; Optype mop = op; - assert(mop != NOP); // silence compiler warning + //assert(mop != NOP); // silence compiler warning if (parent) { merge_point = &merge; @@ -300,13 +523,19 @@ ParanoidNumber * ParanoidNumber::OperationTerm(ParanoidNumber * b, Optype op, Pa for (auto prev : m_next[invop]) { if (prev->OperationTerm(b, rev, merge_point, merge_op) == b) + { + //assert(SanityCheck()); return b; + } } for (auto next : m_next[op]) { if (next->OperationTerm(b, fwd, merge_point, merge_op) == b) + { + //assert(SanityCheck()); return b; + } } @@ -316,6 +545,9 @@ ParanoidNumber * ParanoidNumber::OperationTerm(ParanoidNumber * b, Optype op, Pa { //merge->m_next[*merge_op].push_back(b); m_next[op].push_back(b); + #ifdef PARANOID_SIZE_LIMIT + m_size += b->m_size; + #endif } else { @@ -325,12 +557,36 @@ ParanoidNumber * ParanoidNumber::OperationTerm(ParanoidNumber * b, Optype op, Pa *merge_op = op; } } + + //assert(SanityCheck()); + return NULL; } ParanoidNumber * ParanoidNumber::OperationFactor(ParanoidNumber * b, Optype op, ParanoidNumber ** merge_point, Optype * merge_op) { - m_cached_result = nan(""); + #ifdef PARANOID_CACHE_RESULTS + m_cache_valid = false; + #endif + #ifdef PARANOID_SIZE_LIMIT + if (m_size + b->m_size >= PARANOID_SIZE_LIMIT) + { + this->operator=(this->Digit()); + if (op == MULTIPLY) + m_value *= b->Digit(); + else + m_value /= b->Digit(); + m_size = 1; + #ifdef PARANOID_CACHE_RESULTS + m_cached_result = m_value; + m_cache_valid = true; + #endif + //Debug("Cut off %p", this); + return b; + + } + #endif + if (Floating() && m_value == 0) { return b; @@ -342,15 +598,20 @@ ParanoidNumber * ParanoidNumber::OperationFactor(ParanoidNumber * b, Optype op, for (int i = 0; i < NOP; ++i) { for (auto n : m_next[i]) - delete n; + delete (n); m_next[i].clear(); swap(m_next[i], b->m_next[i]); } + //assert(SanityCheck()); return b; } if (b->Floating() && b->m_value == 1) return b; - + if (b->Floating() && b->m_value == 0 && op == MULTIPLY) + { + operator=(*b); + return b; + } if (NoTerms() && b->NoTerms()) @@ -360,17 +621,14 @@ ParanoidNumber * ParanoidNumber::OperationFactor(ParanoidNumber * b, Optype op, Optype mulop = (op == MULTIPLY) ? MULTIPLY : DIVIDE; for (auto mul : b->m_next[MULTIPLY]) { - delete OperationFactor(mul, mulop); + delete(OperationFactor(mul, mulop)); } Optype divop = (op == MULTIPLY) ? DIVIDE : MULTIPLY; for (auto div : b->m_next[DIVIDE]) - delete OperationFactor(div, divop); + delete(OperationFactor(div, divop)); b->m_next[DIVIDE].clear(); b->m_next[MULTIPLY].clear(); - - - return b; } } @@ -399,23 +657,17 @@ ParanoidNumber * ParanoidNumber::OperationFactor(ParanoidNumber * b, Optype op, rev = DIVIDE; } - ParanoidNumber * cpy_b = NULL; - - if (m_next[ADD].size() > 0 || m_next[SUBTRACT].size() > 0) - { - cpy_b = new ParanoidNumber(*b); - } - + ParanoidNumber * cpy_b = new ParanoidNumber(*b); for (auto prev : m_next[invop]) { if (prev->OperationFactor(b, rev, merge_point, merge_op) == b) { for (auto add : m_next[ADD]) - delete add->OperationFactor(new ParanoidNumber(*cpy_b), op); + delete(add->OperationFactor(new ParanoidNumber(*cpy_b), op)); for (auto sub : m_next[SUBTRACT]) - delete sub->OperationFactor(new ParanoidNumber(*cpy_b), op); + delete(sub->OperationFactor(new ParanoidNumber(*cpy_b), op)); - delete cpy_b; + delete(cpy_b); return b; } } @@ -424,22 +676,35 @@ ParanoidNumber * ParanoidNumber::OperationFactor(ParanoidNumber * b, Optype op, if (next->OperationFactor(b, fwd, merge_point, merge_op) == b) { for (auto add : m_next[ADD]) - delete add->OperationFactor(new ParanoidNumber(*cpy_b), op); + { + delete(add->OperationFactor(new ParanoidNumber(*cpy_b), op)); + } for (auto sub : m_next[SUBTRACT]) - delete sub->OperationFactor(new ParanoidNumber(*cpy_b), op); - delete cpy_b; + { + delete(sub->OperationFactor(new ParanoidNumber(*cpy_b), op)); + } + delete(cpy_b); return b; } } if (parent) { + //assert(b != NULL); m_next[op].push_back(b); for (auto add : m_next[ADD]) - delete add->OperationFactor(new ParanoidNumber(*cpy_b), op); + delete(add->OperationFactor(new ParanoidNumber(*cpy_b), op)); for (auto sub : m_next[SUBTRACT]) - delete sub->OperationFactor(new ParanoidNumber(*cpy_b), op); + delete(sub->OperationFactor(new ParanoidNumber(*cpy_b), op)); + + #ifdef PARANOID_SIZE_LIMIT + m_size += b->m_size; + #endif } + //assert(SanityCheck()); + + + return NULL; } @@ -484,12 +749,85 @@ string ParanoidNumber::PStr() const bool ParanoidNumber::Simplify(Optype op) { - vector next(0); + + if (Floating()) + return false; + + //assert(SanityCheck()); + vector next; + next.clear(); swap(m_next[op], next); + m_next[op].clear(); + //assert(m_next[op].size() == 0); + //assert(SanityCheck()); + Optype fwd = op; + if (op == DIVIDE) + fwd = MULTIPLY; + else if (op == SUBTRACT) + fwd = ADD; + + + vector hold[2]; + if (op == MULTIPLY || op == DIVIDE) + { + swap(m_next[ADD], hold[0]); + swap(m_next[SUBTRACT], hold[1]); + } + + for (vector::iterator n = next.begin(); n != next.end(); ++n) + { + if (*n == NULL) + continue; + for (vector::iterator m = n; m != next.end(); ++m) + { + if ((*m) == (*n)) + continue; + if (*m == NULL) + continue; + + ParanoidNumber * parent = this; + Optype mop = op; + ParanoidNumber * result = (*n)->Operation(*m, fwd, &parent, &mop); + if (result != NULL) + { + #ifdef PARANOID_SIZE_LIMIT + m_size -= result->m_size; + #endif + *m = NULL; + delete(result); + } + } + } + + + for (auto n : next) { - delete Operation(n, op); + if (n != NULL) + { + #ifdef PARANOID_SIZE_LIMIT + if (Operation(n, op) == n) + { + m_size -= n->m_size; + delete n; + } + #else + delete(Operation(n, op)); + #endif + } + } + + if (op == MULTIPLY || op == DIVIDE) + { + swap(m_next[ADD], hold[0]); + swap(m_next[SUBTRACT], hold[1]); } + + set s; + //if (!SanityCheck(s)) + //{ + // Error("Simplify broke Sanity"); + //} return (next.size() > m_next[op].size()); } @@ -503,29 +841,41 @@ bool ParanoidNumber::FullSimplify() return result; } - -ParanoidNumber::digit_t ParanoidNumber::Digit() +ParanoidNumber::digit_t ParanoidNumber::Digit() const { - if (!isnan(m_cached_result)) + + // Get around the absurd requirement that const correctness be observed. + #ifdef PARANOID_CACHE_RESULTS + if (m_cache_valid) // le sigh ambiguous function compiler warnings return m_cached_result; - m_cached_result = m_value; + + digit_t & result = ((ParanoidNumber*)(this))->m_cached_result; + + #else + digit_t result; + #endif + result = m_value; for (auto mul : m_next[MULTIPLY]) { - m_cached_result *= mul->Digit(); + result *= mul->Digit(); } for (auto div : m_next[DIVIDE]) { - m_cached_result /= div->Digit(); + result /= div->Digit(); } for (auto add : m_next[ADD]) - m_cached_result += add->Digit(); + result += add->Digit(); for (auto sub : m_next[SUBTRACT]) - m_cached_result -= sub->Digit(); - return m_cached_result; + result -= sub->Digit(); + + #ifdef PARANOID_CACHE_RESULTS + ((ParanoidNumber*)(this))->m_cache_valid = true; + #endif + return result; } -ParanoidNumber::digit_t ParanoidNumber::GetFactors() +ParanoidNumber::digit_t ParanoidNumber::GetFactors() const { digit_t value = 1; for (auto mul : m_next[MULTIPLY]) @@ -536,7 +886,7 @@ ParanoidNumber::digit_t ParanoidNumber::GetFactors() } -ParanoidNumber::digit_t ParanoidNumber::GetTerms() +ParanoidNumber::digit_t ParanoidNumber::GetTerms() const { digit_t value = 0; for (auto add : m_next[ADD]) @@ -546,5 +896,110 @@ ParanoidNumber::digit_t ParanoidNumber::GetTerms() return value; } +bool ParanoidNumber::SanityCheck(set & visited) const +{ + if (this == NULL) + { + Error("NULL pointer in tree"); + return false; + } + + if (visited.find((ParanoidNumber*)this) != visited.end()) + { + Error("I think I've seen this tree before..."); + return false; + } + + visited.insert((ParanoidNumber*)this); + + for (auto add : m_next[ADD]) + { + if (!add->SanityCheck(visited)) + return false; + } + for (auto sub : m_next[SUBTRACT]) + { + if (!sub->SanityCheck(visited)) + return false; + } + for (auto mul : m_next[MULTIPLY]) + { + if (!mul->SanityCheck(visited)) + return false; + } + + for (auto div : m_next[DIVIDE]) + { + if (!div->SanityCheck(visited)) + return false; + if (div->Digit() == 0) + { + Error("Divide by zero"); + return false; + } + } + return true; +} + +void ParanoidNumber::Negate() +{ + swap(m_next[ADD], m_next[SUBTRACT]); + m_value = -m_value; + #ifdef PARANOID_CACHE_RESULTS + m_cached_result = -m_cached_result; + #endif +} + +#ifdef PARANOID_USE_ARENA + +void * ParanoidNumber::operator new(size_t s) +{ + return g_arena.allocate(s); +} + +void ParanoidNumber::operator delete(void * p) +{ + g_arena.deallocate(p); +} + +ParanoidNumber::Arena::Arena(int64_t block_size) : m_block_size(block_size), m_spare(NULL) +{ + m_blocks.push_back({malloc(block_size*sizeof(ParanoidNumber)),0}); +} + +ParanoidNumber::Arena::~Arena() +{ + for (auto block : m_blocks) + { + free(block.memory); + } + m_blocks.clear(); +} + +void * ParanoidNumber::Arena::allocate(size_t s) +{ + if (m_spare != NULL) + { + void * result = m_spare; + m_spare = NULL; + return result; + } + + Block & b = m_blocks.back(); + void * result = (ParanoidNumber*)(b.memory) + (b.used++); + if (b.used >= m_block_size) + { + m_block_size *= 2; + Debug("Add block of size %d", m_block_size); + m_blocks.push_back({malloc(m_block_size*sizeof(ParanoidNumber)), 0}); + } + return result; +} + +void ParanoidNumber::Arena::deallocate(void * p) +{ + m_spare = p; +} +#endif //PARANOID_USE_ARENA }