diff --git a/barry.hpp b/barry.hpp index d4898cb59..f901730f9 100644 --- a/barry.hpp +++ b/barry.hpp @@ -31,7 +31,7 @@ /* Versioning */ #define BARRY_VERSION_MAYOR 0 #define BARRY_VERSION_MINOR 2 -#define BARRY_VERSION_PATCH 1 +#define BARRY_VERSION_PATCH 2 #define BARRY_VERSION BARRY_VERSION_MAYOR ## . ## BARRY_VERSION_MINOR ## . ## BARRY_VERSION_PATCH static const int barry_version_major = BARRY_VERSION_MAYOR; @@ -396,19 +396,40 @@ struct vecHasher std::size_t operator()(std::vector< T > const& dat) const noexcept { - + + // Well-defined for empty input (the loop below reads dat[0u]). + if (dat.empty()) + return 0u; + std::hash< T > hasher; - std::size_t hash = hasher(dat[0u]); - + + // Finalizer from splitmix64: std::hash of arithmetic types is + // usually the identity (or a bit-cast for doubles, where small + // integers leave the low ~50 bits zero). Without this avalanche + // step the low bits of the combined hash carry almost no entropy, + // and unordered containers that mask the hash with a power-of-two + // bucket count put (nearly) all keys in a single bucket, turning + // lookups into O(n). The mixing runs in an explicit 64-bit + // accumulator so the avalanche is identical where size_t is 32-bit + // (the constants and shifts are 64-bit). + auto mix64 = [](std::uint64_t z) -> std::uint64_t { + z += 0x9e3779b97f4a7c15ULL; + z = (z ^ (z >> 30)) * 0xbf58476d1ce4e5b9ULL; + z = (z ^ (z >> 27)) * 0x94d049bb133111ebULL; + return z ^ (z >> 31); + }; + + std::uint64_t hash = mix64(hasher(dat[0u])); + // ^ makes bitwise XOR // 0x9e3779b9 is a 32 bit constant (comes from the golden ratio) // << is a shift operator, something like lhs * 2^(rhs) if (dat.size() > 1u) for (size_t i = 1u; i < dat.size(); ++i) - hash ^= hasher(dat[i]) + 0x9e3779b9 + (hash<<6) + (hash>>2); - - return hash; - + hash ^= mix64(hasher(dat[i])) + 0x9e3779b9 + (hash<<6) + (hash>>2); + + return static_cast< std::size_t >(hash); + } }; @@ -668,7 +689,7 @@ inline double vec_inner_prod( * - term k * */ -template +template class FreqTable { private: @@ -677,8 +698,56 @@ class FreqTable { size_t k = 0u; size_t n = 0u; + /** + * Rows whose hash collided with a row already in `index` + * (same 64-bit hash, different statistics). Keyed by the shared + * hash, holding the positions of the extra rows in `data`. + * Empty unless a true hash collision occurs, so it adds no + * per-row overhead in the common case. + */ + std::unordered_map> collided; + typename std::unordered_map::iterator iter; - + + /** + * Whether the row stored at `pos` (position of its weight in + * `data`) holds exactly the statistics in `x`. + */ + bool row_matches(size_t pos, const std::vector< T > & x) const + { + for (size_t j = 0u; j < k; ++j) + if (data[pos + 1u + j] != x[j]) + return false; + return true; + } + + /** + * Fallback for a true hash collision: the row at `pos` shares the + * hash `h` with `x` but holds different statistics. Looks `x` up in + * (or appends it to) the per-hash overflow chain. + */ + size_t add_collided(size_t h, const std::vector< T > & x) + { + + auto & chain = collided[h]; + for (auto p : chain) + { + if (row_matches(p, x)) + { + data[p] += 1.0; + return h; + } + } + + chain.push_back(data.size()); + data.push_back(1.0); + data.insert(data.end(), x.begin(), x.end()); + n++; + + return h; + + } + public: // size_t ncols; FreqTable() {}; @@ -739,16 +808,23 @@ inline size_t FreqTable::add( if (x.size() != k) throw std::length_error( - "The value you are trying to add doesn't have the same lenght used in the database." + "The value you are trying to add doesn't have the same length used in the database." ); - + #if __cplusplus > 201700L auto iter2 = index.try_emplace(h, data.size()); - + if (!iter2.second) { - - data[(iter2.first)->second] += 1.0; + + // The hash existed: make sure the row it points to actually + // holds the same statistics before counting it. Otherwise it + // is a hash collision and the row goes to the overflow chain. + size_t pos = (iter2.first)->second; + if (row_matches(pos, x)) + data[pos] += 1.0; + else + return add_collided(h, x); } else @@ -762,25 +838,31 @@ inline size_t FreqTable::add( if (iter == index.end()) { - + index.insert({h, data.size()}); data.push_back(1.0); data.insert(data.end(), x.begin(), x.end()); n++; - + return h; } - data[(*iter).second] += 1.0; - + // The hash existed: make sure the row it points to actually + // holds the same statistics before counting it. Otherwise it + // is a hash collision and the row goes to the overflow chain. + if (row_matches((*iter).second, x)) + data[(*iter).second] += 1.0; + else + return add_collided(h, x); + #endif - + } - + return h; } @@ -819,6 +901,7 @@ inline void FreqTable::clear() { index.clear(); + collided.clear(); data.clear(); n = 0u; @@ -883,7 +966,9 @@ template inline size_t FreqTable::size() const noexcept { - return index.size(); + // Not index.size(): rows stored in the collision overflow chain + // have no entry of their own in -index-. + return n; } @@ -891,17 +976,37 @@ template inline size_t FreqTable::make_hash(const std::vector< T > & x) const { + // Well-defined for empty input (the loop below reads x[0u]). + if (x.empty()) + return 0u; + std::hash< T > hasher; - std::size_t hash = hasher(x[0u]); - + + // Finalizer from splitmix64: std::hash of arithmetic types is usually + // the identity (or a bit-cast for doubles, where small integers leave + // the low ~50 bits zero). Without this avalanche step, the low bits of + // the combined hash carry almost no entropy, and unordered containers + // that mask the hash with a power-of-two bucket count put (nearly) all + // keys in a single bucket, turning lookups into O(n). The mixing runs + // in an explicit 64-bit accumulator so the avalanche is identical where + // size_t is 32-bit (the constants and shifts are 64-bit). + auto mix64 = [](std::uint64_t z) -> std::uint64_t { + z += 0x9e3779b97f4a7c15ULL; + z = (z ^ (z >> 30)) * 0xbf58476d1ce4e5b9ULL; + z = (z ^ (z >> 27)) * 0x94d049bb133111ebULL; + return z ^ (z >> 31); + }; + + std::uint64_t hash = mix64(hasher(x[0u])); + // ^ makes bitwise XOR // 0x9e3779b9 is a 32 bit constant (comes from the golden ratio) // << is a shift operator, something like lhs * 2^(rhs) if (x.size() > 1u) for (size_t i = 1u; i < x.size(); ++i) - hash ^= hasher(x[i]) + 0x9e3779b9 + (hash<<6) + (hash>>2); - - return hash; + hash ^= mix64(hasher(x[i])) + 0x9e3779b9 + (hash<<6) + (hash>>2); + + return static_cast< size_t >(hash); } diff --git a/defm.hpp b/defm.hpp index 361486c8b..5550a9d40 100644 --- a/defm.hpp +++ b/defm.hpp @@ -286,9 +286,24 @@ typedef barry::Rules DEFMRulesDyn; * * ## Transition effects * - * Transition effects can be specified using two sets of curly brackets and - * an greater-than symbol, i.e., `{...} > {...}`. The first set of brackets, - * which we call LHS, can only hold `row id` that are less than `m_order`. + * Transition effects can be specified using curly brackets separated by + * greater-than symbols ('>'). + * + * **Two-group mode (backwards compatible):** `{...} > {...}` + * - First group (LHS): Variables at times 0 to m_order-1. When m_order > 1, + * row indices must be explicitly specified. + * - Second group (RHS): Variables at time m_order. Row indices can be omitted + * and will default to m_order. + * + * **Multi-group mode (explicit):** `{...} > {...} > ... > {...}` (m_order+1 groups) + * - Each group corresponds to a specific time point (0, 1, ..., m_order). + * - Row indices can be omitted and will be inferred from group position. + * - If specified, row indices must match the group position. + * + * Examples: + * - Order 1: `{y0_0} > {y0_1}` or `{y0_0} > {y0}` (both valid) + * - Order 2: `{y0_0} > {y0}` (2-group, implicit final time) + * - Order 2: `{y0_0} > {y0_1} > {y0_2}` (3-group, all explicit) * * * @param formula A string specifying the motif formula (see details). @@ -318,9 +333,10 @@ inline void defm_motif_parser( std::string("\\{\\s*[01]?y[0-9]+(_[0-9]+)?(\\s*,\\s*[01]?y[0-9]+(_[0-9]+)?)*\\s*\\}") + std::string("(\\s*x\\s*[^\\s]+([(].+[)])?\\s*)?") ); + // Updated pattern to match one or more bracketed groups separated by '>' std::regex pattern_transition( - std::string("\\{\\s*[01]?y[0-9]+(_[0-9]+)?(\\s*,\\s*[01]?y[0-9]+(_[0-9]+)?)*\\}\\s*(>)\\s*") + std::string("\\{\\s*[01]?y[0-9]+(_[0-9]+)?(\\s*,\\s*[01]?y[0-9]+(_[0-9]+)?)*\\s*\\}") + + std::string("(\\s*>\\s*\\{\\s*[01]?y[0-9]+(_[0-9]+)?(\\s*,\\s*[01]?y[0-9]+(_[0-9]+)?)*\\s*\\})+") + std::string("(\\s*x\\s*[^\\s]+([(].+[)])?\\s*)?") ); @@ -358,79 +374,187 @@ inline void defm_motif_parser( } - // Will indicate where the arrow is located at - size_t arrow_position = match.position(4u); - - // This pattern will match - std::regex pattern("(0?)y([0-9]+)(_([0-9]+))?"); + // Find all bracketed groups to determine which time point each variable belongs to + std::regex bracket_pattern("\\{[^}]+\\}"); + std::vector> bracket_ranges; // start, end positions + + auto brackets_begin = std::sregex_iterator(formula.begin(), formula.end(), bracket_pattern); + for (auto i = brackets_begin; i != empty; ++i) + bracket_ranges.push_back({i->position(), i->position() + i->length()}); - auto iter = std::sregex_iterator(formula.begin(), formula.end(), pattern); + size_t num_groups = bracket_ranges.size(); - for (auto i = iter; i != empty; ++i) + // For backwards compatibility, allow 2 groups (original behavior) or m_order+1 groups (new behavior) + if (num_groups == 2) { + // Two-group mode (backwards compatible): + // - First group: variables at time 0 to m_order-1 (must have explicit row when m_order > 1) + // - Second group: variables at time m_order (row can be implicit or explicit) - // Baseline position - size_t current_location = i->position(0u); - - // First value true/false - bool is_positive = true; - if (i->operator[](1u).str() == "0") - is_positive = false; + // This pattern will match + std::regex pattern("(0?)y([0-9]+)(_([0-9]+))?"); - // Variable position - size_t y_col = std::stoul(i->operator[](2u).str()); - if (y_col >= y_ncol) - throw std::logic_error("The proposed column is out of range."); + auto iter = std::sregex_iterator(formula.begin(), formula.end(), pattern); - // Time location - size_t y_row; - std::string tmp_str = i->operator[](4u).str(); - if (m_order > 1) + for (auto i = iter; i != empty; ++i) { - // If missing, we replace with the location - if (tmp_str == "") + + // Baseline position + size_t current_location = i->position(0u); + + // First value true/false + bool is_positive = true; + if (i->operator[](1u).str() == "0") + is_positive = false; + + // Variable position + size_t y_col = std::stoul(i->operator[](2u).str()); + if (y_col >= y_ncol) + throw std::logic_error("The proposed column is out of range."); + + // Time location + size_t y_row; + std::string tmp_str = i->operator[](4u).str(); + if (m_order > 1) { + // If missing, we replace with the location + if (tmp_str == "") + { + + if (current_location >= bracket_ranges[1].first) + y_row = m_order; + else + throw std::logic_error("LHS of transition must specify time when m_order > 1"); + + } else + y_row = std::stoul(tmp_str); + + if (y_row > m_order) + throw std::logic_error("The proposed row is out of range."); - if (current_location > arrow_position) - y_row = m_order; + + } else { + + // If missing, we replace with the location + if (tmp_str != "") + y_row = std::stoul(tmp_str); else - throw std::logic_error("LHS of transition must specify time when m_order > 1"); + y_row = (current_location < bracket_ranges[1].first ? 0u: 1u); - } else - y_row = std::stoul(tmp_str); + } - if (y_row > m_order) - throw std::logic_error("The proposed row is out of range."); + if (selected[y_col * (m_order + 1) + y_row]) + throw std::logic_error( + "The term " + i->str() + " shows more than once in the formula."); + + // Only variables at time m_order can be in the RHS (second bracketed group) + if ((current_location >= bracket_ranges[1].first) && (y_row != m_order)) + throw std::logic_error( + "Only the row " + std::to_string(m_order) + + " can be specified at the RHS of the motif." + ); + + selected[y_col * (m_order + 1) + y_row] = true; + + locations.push_back(y_col * (m_order + 1) + y_row); + signs.push_back(is_positive); + + + } + } + else if (num_groups == m_order + 1) + { + // New behavior: each group corresponds to a time point (0 to m_order) + + // This pattern will match + std::regex pattern("(0?)y([0-9]+)(_([0-9]+))?"); + auto iter = std::sregex_iterator(formula.begin(), formula.end(), pattern); - } else { + for (auto i = iter; i != empty; ++i) + { - // If missing, we replace with the location + // Baseline position + size_t current_location = i->position(0u); + + // First value true/false + bool is_positive = true; + if (i->operator[](1u).str() == "0") + is_positive = false; + + // Variable position + size_t y_col = std::stoul(i->operator[](2u).str()); + if (y_col >= y_ncol) + throw std::logic_error("The proposed column is out of range."); + + // Determine which bracketed group this variable belongs to + // Note: The regex pattern ensures all variables are within bracketed groups, + // as the pattern only matches content inside brackets. Variables in covariate + // expressions (after 'x') are handled separately by pattern_conditional above. + size_t group_idx = 0; + bool found_group = false; + for (size_t g = 0; g < bracket_ranges.size(); ++g) + { + if (current_location >= bracket_ranges[g].first && + current_location < bracket_ranges[g].second) + { + group_idx = g; + found_group = true; + break; + } + } + + // Safety check: ensure the variable was found in a bracketed group + // This should never happen given the regex, but verify for safety + if (!found_group) + throw std::logic_error( + "Internal error: variable " + i->str() + + " not found within any bracketed group."); + + // Time location + size_t y_row; + std::string tmp_str = i->operator[](4u).str(); + + // If row is explicitly specified, use it; otherwise infer from group position if (tmp_str != "") + { y_row = std::stoul(tmp_str); + + // Validate that explicit row matches the expected group position + if (y_row != group_idx) + throw std::logic_error( + "Explicit row index " + std::to_string(y_row) + + " does not match the position in the formula (group " + + std::to_string(group_idx) + ")."); + } else - y_row = (current_location < arrow_position ? 0u: 1u); - - } + { + // Infer row from bracketed group position + y_row = group_idx; + } - if (selected[y_col * (m_order + 1) + y_row]) - throw std::logic_error( - "The term " + i->str() + " shows more than once in the formula."); + if (y_row > m_order) + throw std::logic_error("The proposed row is out of range."); - // Only the end of the chain can be located at position after the - // arrow - if ((current_location > arrow_position) && (y_row != m_order)) - throw std::logic_error( - "Only the row " + std::to_string(m_order) + - " can be specified at the RHS of the motif." - ); + if (selected[y_col * (m_order + 1) + y_row]) + throw std::logic_error( + "The term " + i->str() + " shows more than once in the formula."); - selected[y_col * (m_order + 1) + y_row] = true; + selected[y_col * (m_order + 1) + y_row] = true; - locations.push_back(y_col * (m_order + 1) + y_row); - signs.push_back(is_positive); - + locations.push_back(y_col * (m_order + 1) + y_row); + signs.push_back(is_positive); + + } + } + else + { + throw std::logic_error( + "For a Markov model of order " + std::to_string(m_order) + + ", transition formulas must have either 2 bracketed groups " + + "(for backwards compatibility) or exactly " + std::to_string(m_order + 1) + + " bracketed groups (found " + std::to_string(num_groups) + ")."); } return; @@ -1007,43 +1131,41 @@ inline void counter_generic( } - // Checking if any prior to the event - bool any_before_event = false; - - for (size_t i = 0u; i < m_order; ++i) + #ifdef BARRY_WITH_LATEX + name += "$"; + #endif + + // Generate name by showing each time point separately + // Loop through all time points (0 to m_order) + for (size_t i = 0u; i <= m_order; ++i) { + // Check if this time point has any variables + bool has_vars = false; for (size_t j = 0u; j < n_y; ++j) { - if (motif(i,j) != 0) + if (motif(i, j) != 0) { - any_before_event = true; + has_vars = true; break; } - } - } - - #ifdef BARRY_WITH_LATEX - name += "$"; - #endif - if (any_before_event) + if (!has_vars) + continue; + + // Open bracket for this time point #ifdef BARRY_WITH_LATEX name += "("; #else name += "{"; #endif - // If order is greater than zero, the starting point of the transtion - for (size_t i = 0u; i < m_order; ++i) - { - + // Add variables for this time point bool row_start = true; for (size_t j = 0u; j < n_y; ++j) { - // Is it included? - if (motif(i,j) == 0) + if (motif(i, j) == 0) continue; // Is not the first? @@ -1058,60 +1180,28 @@ inline void counter_generic( name += (std::string("y") + std::to_string(j)); #ifdef BARRY_WITH_LATEX - name += (motif(i,j) < 0 ? "^-" : "^+"); + name += (motif(i, j) < 0 ? "^-" : "^+"); #else - name += (motif(i,j) < 0 ? "-" : "+"); + name += (motif(i, j) < 0 ? "-" : "+"); #endif - } - } - - // If it has starting point, then need to close. - if (any_before_event & (m_order > 0u)) - #ifdef BARRY_WITH_LATEX - name += ") -> ("; - #else - name += std::string("}") + u8"\u21E8" + std::string("{"); - #endif - else + // Close bracket for this time point #ifdef BARRY_WITH_LATEX - name += "("; + name += ")"; #else - name += "{"; + name += "}"; #endif - // Looking onto the transtions - bool row_start = true; - for (size_t j = 0u; j < n_y; ++j) - { - - if (motif(m_order, j) == 0) - continue; - - if (row_start) - row_start = false; - else - name += ", "; - - if (y_names != nullptr) - name += y_names->operator[](j); - else - name += (std::string("y") + std::to_string(j)); - - #ifdef BARRY_WITH_LATEX - name += (motif(m_order, j) < 0 ? "^-" : "^+" ); - #else - name += (motif(m_order, j) < 0 ? "-" : "+" ); - #endif - - + // Add arrow if not the last time point + if (i < m_order) + { + name += ">"; + } } #ifdef BARRY_WITH_LATEX - name += ")$"; - #else - name += "}"; + name += "$"; #endif if (covar_index >= 0) @@ -1345,7 +1435,7 @@ inline void rule_constrain_support( support->get_counters()->get_names()[term_index] + "' within [" + std::to_string(lb) + ", " + std::to_string(ub) + std::string("]"), - std::string("When the support is ennumerated, only states where the statistic '") + + std::string("When the support is enumerated, only states where the statistic '") + support->get_counters()->get_names()[term_index] + std::to_string(term_index) + "' falls within [" + std::to_string(lb) + ", " + std::to_string(ub) + "] are included." diff --git a/include/barry/barry.hpp b/include/barry/barry.hpp index c76a44935..a80aa6a86 100644 --- a/include/barry/barry.hpp +++ b/include/barry/barry.hpp @@ -31,7 +31,7 @@ /* Versioning */ #define BARRY_VERSION_MAYOR 0 #define BARRY_VERSION_MINOR 2 -#define BARRY_VERSION_PATCH 1 +#define BARRY_VERSION_PATCH 2 #define BARRY_VERSION BARRY_VERSION_MAYOR ## . ## BARRY_VERSION_MINOR ## . ## BARRY_VERSION_PATCH static const int barry_version_major = BARRY_VERSION_MAYOR; diff --git a/include/barry/freqtable.hpp b/include/barry/freqtable.hpp index 269978e7f..d958eba88 100644 --- a/include/barry/freqtable.hpp +++ b/include/barry/freqtable.hpp @@ -18,7 +18,7 @@ * - term k * */ -template +template class FreqTable { private: @@ -27,8 +27,56 @@ class FreqTable { size_t k = 0u; size_t n = 0u; + /** + * Rows whose hash collided with a row already in `index` + * (same 64-bit hash, different statistics). Keyed by the shared + * hash, holding the positions of the extra rows in `data`. + * Empty unless a true hash collision occurs, so it adds no + * per-row overhead in the common case. + */ + std::unordered_map> collided; + typename std::unordered_map::iterator iter; - + + /** + * Whether the row stored at `pos` (position of its weight in + * `data`) holds exactly the statistics in `x`. + */ + bool row_matches(size_t pos, const std::vector< T > & x) const + { + for (size_t j = 0u; j < k; ++j) + if (data[pos + 1u + j] != x[j]) + return false; + return true; + } + + /** + * Fallback for a true hash collision: the row at `pos` shares the + * hash `h` with `x` but holds different statistics. Looks `x` up in + * (or appends it to) the per-hash overflow chain. + */ + size_t add_collided(size_t h, const std::vector< T > & x) + { + + auto & chain = collided[h]; + for (auto p : chain) + { + if (row_matches(p, x)) + { + data[p] += 1.0; + return h; + } + } + + chain.push_back(data.size()); + data.push_back(1.0); + data.insert(data.end(), x.begin(), x.end()); + n++; + + return h; + + } + public: // size_t ncols; FreqTable() {}; @@ -89,16 +137,23 @@ inline size_t FreqTable::add( if (x.size() != k) throw std::length_error( - "The value you are trying to add doesn't have the same lenght used in the database." + "The value you are trying to add doesn't have the same length used in the database." ); - + #if __cplusplus > 201700L auto iter2 = index.try_emplace(h, data.size()); - + if (!iter2.second) { - - data[(iter2.first)->second] += 1.0; + + // The hash existed: make sure the row it points to actually + // holds the same statistics before counting it. Otherwise it + // is a hash collision and the row goes to the overflow chain. + size_t pos = (iter2.first)->second; + if (row_matches(pos, x)) + data[pos] += 1.0; + else + return add_collided(h, x); } else @@ -112,25 +167,31 @@ inline size_t FreqTable::add( if (iter == index.end()) { - + index.insert({h, data.size()}); data.push_back(1.0); data.insert(data.end(), x.begin(), x.end()); n++; - + return h; } - data[(*iter).second] += 1.0; - + // The hash existed: make sure the row it points to actually + // holds the same statistics before counting it. Otherwise it + // is a hash collision and the row goes to the overflow chain. + if (row_matches((*iter).second, x)) + data[(*iter).second] += 1.0; + else + return add_collided(h, x); + #endif - + } - + return h; } @@ -169,6 +230,7 @@ inline void FreqTable::clear() { index.clear(); + collided.clear(); data.clear(); n = 0u; @@ -233,7 +295,9 @@ template inline size_t FreqTable::size() const noexcept { - return index.size(); + // Not index.size(): rows stored in the collision overflow chain + // have no entry of their own in -index-. + return n; } @@ -241,17 +305,37 @@ template inline size_t FreqTable::make_hash(const std::vector< T > & x) const { + // Well-defined for empty input (the loop below reads x[0u]). + if (x.empty()) + return 0u; + std::hash< T > hasher; - std::size_t hash = hasher(x[0u]); - + + // Finalizer from splitmix64: std::hash of arithmetic types is usually + // the identity (or a bit-cast for doubles, where small integers leave + // the low ~50 bits zero). Without this avalanche step, the low bits of + // the combined hash carry almost no entropy, and unordered containers + // that mask the hash with a power-of-two bucket count put (nearly) all + // keys in a single bucket, turning lookups into O(n). The mixing runs + // in an explicit 64-bit accumulator so the avalanche is identical where + // size_t is 32-bit (the constants and shifts are 64-bit). + auto mix64 = [](std::uint64_t z) -> std::uint64_t { + z += 0x9e3779b97f4a7c15ULL; + z = (z ^ (z >> 30)) * 0xbf58476d1ce4e5b9ULL; + z = (z ^ (z >> 27)) * 0x94d049bb133111ebULL; + return z ^ (z >> 31); + }; + + std::uint64_t hash = mix64(hasher(x[0u])); + // ^ makes bitwise XOR // 0x9e3779b9 is a 32 bit constant (comes from the golden ratio) // << is a shift operator, something like lhs * 2^(rhs) if (x.size() > 1u) for (size_t i = 1u; i < x.size(); ++i) - hash ^= hasher(x[i]) + 0x9e3779b9 + (hash<<6) + (hash>>2); - - return hash; + hash ^= mix64(hasher(x[i])) + 0x9e3779b9 + (hash<<6) + (hash>>2); + + return static_cast< size_t >(hash); } diff --git a/include/barry/typedefs.hpp b/include/barry/typedefs.hpp index 5d34c90bf..c59670528 100644 --- a/include/barry/typedefs.hpp +++ b/include/barry/typedefs.hpp @@ -107,19 +107,40 @@ struct vecHasher std::size_t operator()(std::vector< T > const& dat) const noexcept { - + + // Well-defined for empty input (the loop below reads dat[0u]). + if (dat.empty()) + return 0u; + std::hash< T > hasher; - std::size_t hash = hasher(dat[0u]); - + + // Finalizer from splitmix64: std::hash of arithmetic types is + // usually the identity (or a bit-cast for doubles, where small + // integers leave the low ~50 bits zero). Without this avalanche + // step the low bits of the combined hash carry almost no entropy, + // and unordered containers that mask the hash with a power-of-two + // bucket count put (nearly) all keys in a single bucket, turning + // lookups into O(n). The mixing runs in an explicit 64-bit + // accumulator so the avalanche is identical where size_t is 32-bit + // (the constants and shifts are 64-bit). + auto mix64 = [](std::uint64_t z) -> std::uint64_t { + z += 0x9e3779b97f4a7c15ULL; + z = (z ^ (z >> 30)) * 0xbf58476d1ce4e5b9ULL; + z = (z ^ (z >> 27)) * 0x94d049bb133111ebULL; + return z ^ (z >> 31); + }; + + std::uint64_t hash = mix64(hasher(dat[0u])); + // ^ makes bitwise XOR // 0x9e3779b9 is a 32 bit constant (comes from the golden ratio) // << is a shift operator, something like lhs * 2^(rhs) if (dat.size() > 1u) for (size_t i = 1u; i < dat.size(); ++i) - hash ^= hasher(dat[i]) + 0x9e3779b9 + (hash<<6) + (hash>>2); - - return hash; - + hash ^= mix64(hasher(dat[i])) + 0x9e3779b9 + (hash<<6) + (hash>>2); + + return static_cast< std::size_t >(hash); + } }; diff --git a/tests/19-freqtable-collisions.cpp b/tests/19-freqtable-collisions.cpp new file mode 100644 index 000000000..1232b48e8 --- /dev/null +++ b/tests/19-freqtable-collisions.cpp @@ -0,0 +1,59 @@ +#include "tests.hpp" + +BARRY_TEST_CASE("FreqTable keeps distinct stats separate on hash collisions", "[freqtable]") { + + // Regular usage: equal vectors merge, distinct vectors get their own row + barry::FreqTable tab; + std::vector< double > a = {1.0, 2.0, 3.0}; + std::vector< double > b = {4.0, 5.0, 6.0}; + + tab.add(a, nullptr); + tab.add(a, nullptr); + tab.add(b, nullptr); + + REQUIRE(tab.size() == 2u); + + auto counts = tab.as_vector(); + REQUIRE(counts[0u].first == a); + REQUIRE(counts[0u].second == 2u); + REQUIRE(counts[1u].first == b); + REQUIRE(counts[1u].second == 1u); + + // Forced collisions: the precomputed-hash argument lets us feed + // different stat vectors under the same 64-bit hash, mimicking a + // genuine make_hash collision. They must remain separate rows. + barry::FreqTable tab2; + std::vector< double > c = {7.0, 8.0, 9.0}; + size_t h = 42u; + + tab2.add(a, &h); + tab2.add(b, &h); + tab2.add(a, &h); + tab2.add(c, &h); + tab2.add(b, &h); + + REQUIRE(tab2.size() == 3u); + + auto counts2 = tab2.as_vector(); + REQUIRE(counts2[0u].first == a); + REQUIRE(counts2[0u].second == 2u); + REQUIRE(counts2[1u].first == b); + REQUIRE(counts2[1u].second == 2u); + REQUIRE(counts2[2u].first == c); + REQUIRE(counts2[2u].second == 1u); + + // Clearing must also reset the collision chains + tab2.clear(); + REQUIRE(tab2.size() == 0u); + + tab2.add(a, &h); + tab2.add(b, &h); + REQUIRE(tab2.size() == 2u); + + // make_hash must be well-defined (no out-of-bounds read) for an empty + // vector; it returns a fixed value rather than reading element 0. + barry::FreqTable tab3; + std::vector< double > empty; + REQUIRE(tab3.make_hash(empty) == tab3.make_hash(empty)); + +} diff --git a/tests/main.cpp b/tests/main.cpp index fb9a3064a..d2b4b3ffa 100644 --- a/tests/main.cpp +++ b/tests/main.cpp @@ -29,3 +29,4 @@ #include "16b-defm-counts-with-formulas.cpp" #include "17-defm-likelihood.cpp" #include "18-defm-counter-names.cpp" +#include "19-freqtable-collisions.cpp"