From 1d764d29da71a11857c22220c75e56d7e0b0693d Mon Sep 17 00:00:00 2001 From: Lance Fletcher Date: Thu, 5 Sep 2024 11:42:01 -0700 Subject: [PATCH 01/15] Added alias table to ygm. Created new namespace ygm::random. --- include/ygm/container/bag.hpp | 6 +- include/ygm/random/alias_table.hpp | 296 ++++++++++++++++++ include/ygm/random/multi_alias_table.hpp | 364 +++++++++++++++++++++++ include/ygm/{ => random}/random.hpp | 2 +- test/CMakeLists.txt | 1 + test/test_alias_table.cpp | 80 +++++ test/test_bag.cpp | 10 +- test/test_gather_topk.cpp | 1 - test/test_random.cpp | 4 +- test/test_reduce.cpp | 1 - test/test_transform.cpp | 1 - 11 files changed, 752 insertions(+), 14 deletions(-) create mode 100644 include/ygm/random/alias_table.hpp create mode 100644 include/ygm/random/multi_alias_table.hpp rename include/ygm/{ => random}/random.hpp (95%) create mode 100644 test/test_alias_table.cpp diff --git a/include/ygm/container/bag.hpp b/include/ygm/container/bag.hpp index 86511f5b..773c7e2e 100644 --- a/include/ygm/container/bag.hpp +++ b/include/ygm/container/bag.hpp @@ -13,7 +13,7 @@ #include #include #include -#include +#include namespace ygm::container { @@ -210,7 +210,7 @@ class bag : public detail::base_async_insert_value, std::tuple>, } void local_shuffle() { - ygm::default_random_engine<> r(m_comm, std::random_device()()); + ygm::random::default_random_engine<> r(m_comm, std::random_device()()); local_shuffle(r); } @@ -231,7 +231,7 @@ class bag : public detail::base_async_insert_value, std::tuple>, } void global_shuffle() { - ygm::default_random_engine<> r(m_comm, std::random_device()()); + ygm::random::default_random_engine<> r(m_comm, std::random_device()()); global_shuffle(r); } diff --git a/include/ygm/random/alias_table.hpp b/include/ygm/random/alias_table.hpp new file mode 100644 index 00000000..4818ccf3 --- /dev/null +++ b/include/ygm/random/alias_table.hpp @@ -0,0 +1,296 @@ +#include +#include +#include +#include +#include +#include +#include +#include +#include + + +namespace ygm::random { + + + +template +concept convertable_to_double = std::convertible_to; + +template +concept is_pair = requires { + typename T::first_type; + typename T::second_type; +} && std::same_as>; + +template +concept is_associative_container = requires { + typename T::key_type; + typename T::mapped_type; +}; + +template +concept is_value_container = requires { + typename T::value_type; +}; + + +template +class alias_table { + + using self_type = alias_table; + using weight_type = double; + + public: + + struct item { + item_type id; + weight_type weight; + + template + void serialize(Archive& ar) { // Needed for cereal serialization + ar(id, weight); + } + }; + + struct table_item { + double p; // prob item a is selected = p, prob item b is selected = (1 - p) + item_type a; + item_type b; + + template + void serialize(Archive& ar) { // Needed for cereal serialization + ar(p, a, b); + } + }; + + + template + requires is_value_container && + is_pair && + convertable_to_double + alias_table(ygm::comm &comm, RNG &rng, Container &c) : + m_comm(comm), + pthis(this), + m_rng(rng), + m_rank_dist(0, comm.size()-1) { + c.for_all([&](const auto& id_weight_pair){ + m_local_items.push_back({id_weight_pair.first, id_weight_pair.second}); + }); + comm.barrier(); + balance_weight(); + comm.barrier(); + YGM_ASSERT_RELEASE(check_balancing()); + build_alias_table(); + m_local_items.clear(); + } + + template + requires is_associative_container && + convertable_to_double + alias_table(ygm::comm &comm, RNG &rng, Container &c) : + m_comm(comm), + pthis(this), + m_rng(rng), + m_rank_dist(0, comm.size()-1) { + c.for_all([&](const auto &key, const auto &value){ + m_local_items.push_back({key, value}); + }); + comm.barrier(); + balance_weight(); + comm.barrier(); + YGM_ASSERT_RELEASE(check_balancing()); + build_alias_table(); + m_local_items.clear(); + } + + // For debugging + std::string get_weight_str(std::vector items) { + std::string weights_str = "< "; + for (auto& itm : items) { + weights_str += std::to_string(itm.weight) + " "; + } + return weights_str + ">"; + } + + void balance_weight() { + + MPI_Comm comm = m_comm.get_mpi_comm(); + double local_weight = 0.0; + for (int i = 0; i < m_local_items.size(); i++) { + local_weight += m_local_items[i].weight; + } + + double global_weight; + MPI_Allreduce(&local_weight, &global_weight, 1, MPI_DOUBLE, MPI_SUM, comm); + + double prfx_sum_weight; + MPI_Exscan(&local_weight, &prfx_sum_weight, 1, MPI_DOUBLE, MPI_SUM, comm); + + + double target_weight = global_weight / m_comm.size(); // target weight per rank + uint32_t dest_rank = prfx_sum_weight / target_weight; + double curr_weight = std::fmod(prfx_sum_weight, target_weight); // Spillage weight + + std::vector new_local_items; + using ygm_items_ptr = ygm::ygm_ptr>; + ygm_items_ptr ptr_new_items = m_comm.make_ygm_ptr(new_local_items); + m_comm.barrier(); // ensure ygm_ptr is created on every rank + + std::vector items_to_send; + for (uint64_t i = 0; i < m_local_items.size(); i++) { // WARNING: size of m_local_items can grow during loop + item local_item = m_local_items[i]; + if (curr_weight + local_item.weight >= target_weight) { + + double remaining_weight = curr_weight + local_item.weight - target_weight; + double weight_to_send = local_item.weight - remaining_weight; + curr_weight += weight_to_send; + item item_to_send = {local_item.id, weight_to_send}; + items_to_send.push_back(item_to_send); + + if ((curr_weight > 0.0001) && (dest_rank < m_comm.size())) { // Accounts for rounding errors + // Moves weights to dest_rank's new_local_items + m_comm.async(dest_rank, [](std::vector items, ygm_items_ptr new_items_ptr) { + new_items_ptr->insert(new_items_ptr->end(), items.begin(), items.end()); + }, items_to_send, ptr_new_items); + } + // Handle case where item weight is large enough to span multiple rank's alias tables + if (remaining_weight >= target_weight) { + m_local_items.push_back({local_item.id, remaining_weight}); + curr_weight = 0; + } else { + curr_weight = remaining_weight; + } + items_to_send.clear(); + if (curr_weight != 0) { + items_to_send.push_back({local_item.id, curr_weight}); + } + dest_rank++; + + } else { + items_to_send.push_back(local_item); + curr_weight += local_item.weight; + } + } + + // Need to handle items left in items to send + if (items_to_send.size() > 0 && curr_weight > 0.00001 && dest_rank < m_comm.size()) { // Account for floating point errors + m_comm.async(dest_rank, [](std::vector items, ygm_items_ptr new_items_ptr) { + new_items_ptr->insert(new_items_ptr->end(), items.begin(), items.end()); + }, items_to_send, ptr_new_items); + } + + m_comm.barrier(); + std::swap(new_local_items, m_local_items); + } + + bool check_balancing() { + // Not a perfect check. Will check if every rank has the same weight for each table + // but does not confirm that the actual target weight is met on each rank + double local_weight = 0.0; + for (const auto& itm : m_local_items) { + local_weight += itm.weight; + } + m_comm.barrier(); + auto equal = [](double a, double b){ + return (std::abs(a - b) < 0.00001); + }; + bool balanced = ygm::is_same(local_weight, m_comm, equal); + return balanced; + } + + void build_alias_table() { + // Make average weight of items 1, so coin flip is simple across all tables + double local_weight = 0.0; + for (const auto& itm : m_local_items) { + local_weight += itm.weight; + } + double avg_weight = local_weight / m_local_items.size(); + + double new_total_weight = 0; + for (auto& itm : m_local_items) { + itm.weight /= avg_weight; + new_total_weight += itm.weight; + } + double new_avg_weight = new_total_weight / m_local_items.size(); // Should be 1 + if (std::abs((new_avg_weight - 1.0)) < 0.0001) { + YGM_ASSERT_RELEASE(std::abs((new_avg_weight - 1.0)) < 0.00001); + } + + // Implementation of Vose's algorithm, utilized Keith Schwarz numerically stable version + // https://www.keithschwarz.com/darts-dice-coins/ + std::vector heavy_items; + std::vector light_items; + for (auto& itm : m_local_items) { + if (itm.weight < 1.0) { + light_items.push_back(itm); + } else { + heavy_items.push_back(itm); + } + } + + while (!light_items.empty() && !heavy_items.empty()) { + item l = light_items.back(); + item h = heavy_items.back(); + table_item tbl_itm = {l.weight, l.id, h.id}; + m_local_alias_table.push_back(tbl_itm); + h.weight = (h.weight + l.weight) - 1; + light_items.pop_back(); + if (h.weight < 1) { + light_items.push_back(h); + heavy_items.pop_back(); + } + } + + // Either heavy items or light_items is empty, need to flush the non empty + // vector and add them to the alias table with a p value of 1 + while (!heavy_items.empty()) { + item h = heavy_items.back(); + table_item tbl_itm = {1.0, h.id, 0}; + m_local_alias_table.push_back(tbl_itm); + heavy_items.pop_back(); + } + while (!light_items.empty()) { + item l = light_items.back(); + table_item tbl_itm = {1.0, l.id, 0}; + m_local_alias_table.push_back(tbl_itm); + light_items.pop_back(); + } + m_comm.barrier(); + } + + template + void async_sample(Visitor visitor, const VisitorArgs &...args) { // Sample should be provided a lambda/functor that takes as an argument the item type + + auto sample_wrapper = [](auto ptr_MAT, const VisitorArgs &...args) { + std::uniform_int_distribution table_item_dist(0, ptr_MAT->m_local_alias_table.size()-1); + table_item tbl_itm = ptr_MAT->m_local_alias_table[table_item_dist(ptr_MAT->m_rng)]; + item_type s; + if (tbl_itm.p == 1) { + s = tbl_itm.a; + } else { + std::uniform_real_distribution zero_one_dist(0.0, 1.0); + float f = zero_one_dist(ptr_MAT->m_rng); + if (f < tbl_itm.p) { + s = tbl_itm.a; + } else { + s = tbl_itm.b; + } + } + Visitor *vis = nullptr; + ygm::meta::apply_optional(*vis, std::make_tuple(ptr_MAT), std::forward_as_tuple(s, args...)); + }; + + uint32_t dest_rank = m_rank_dist(m_rng); + m_comm.async(dest_rank, sample_wrapper, pthis, std::forward(args)...); + } + + private: + ygm::comm& m_comm; + ygm::ygm_ptr pthis; + std::uniform_int_distribution m_rank_dist; + RNG& m_rng; + std::vector m_local_items; + std::vector m_local_alias_table; +}; + +} // namespace ygm::random \ No newline at end of file diff --git a/include/ygm/random/multi_alias_table.hpp b/include/ygm/random/multi_alias_table.hpp new file mode 100644 index 00000000..4faeb92d --- /dev/null +++ b/include/ygm/random/multi_alias_table.hpp @@ -0,0 +1,364 @@ +#include +#include +#include +#include +#include +#include +#include +#include +#include + + +namespace ygm::random { + +template +class multi_alias_table { + + using self_type = multi_alias_table; + + public: + + struct item { + table_id_type t_id; + item_id_type i_id; + weight_type weight; + + template + void serialize(Archive& ar) { // Needed for cereal serialization + ar(t_id, i_id, weight); + } + }; + + + struct table_item { + double p; // prob item a is selected = p, prob item b is selected = (1 - p) + item_id_type a; + item_id_type b; + + template + void serialize(Archive& ar) { // Needed for cereal serialization + ar(p, a, b); + } + }; + + + multi_alias_table(ygm::comm &comm, RNG rng) : + m_comm(comm), + pthis(this), + m_rng(rng), + m_tables_built(false), + m_balanced(false), + // m_local_weight(0), + m_rank_dist(0, comm.size()-1) { } + + + template + void load_items_from_container(container &cont) { + cont.for_all([&](item itm){ + m_local_items[itm.t_id].push_back(itm); + m_each_tables_local_weight[itm.t_id] += itm.weight; + }); + } + + void local_add_item (std::tuple i) { + item itm({std::get<0>(i), std::get<1>(i), std::get<2>(i)}); + m_local_items[itm.t_id].push_back(itm); + m_each_tables_local_weight[itm.t_id] += itm.weight; + } + + void add_item(std::tuple i) { + item itm({std::get<0>(i), std::get<1>(i), std::get<2>(i)}); + + auto add_item_to_rand_rank = [](item i, auto ptr_MAT){ + ptr_MAT->m_local_items[i.t_id].push_back(i); + ptr_MAT->m_each_tables_local_weight[i.t_id] += i.weight; + }; + + m_comm.async(m_rank_dist(m_rng), add_item_to_rand_rank, itm, pthis); + } + + // For debugging + std::string get_weight_str(std::vector items) { + std::string weights_str = "< "; + for (auto& itm : items) { + weights_str += std::to_string(itm.weight) + " "; + } + return weights_str + ">"; + } + + std::vector get_table_ids() { + std::vector table_ids; + for (const table_id_type t : m_each_tables_local_weight) { + table_ids.push_back(t); + } + return table_ids; + } + + void clear() { + m_comm.barrier(); + m_tables_built = false; + m_local_items.clear(); + m_local_alias_tables.clear(); + m_each_tables_local_weight.clear(); + } + + void clear_items() { + m_local_items.clear(); + } + + void balance_weight() { + using t_id_vec = std::vector; + t_id_vec table_ids; + for (const auto& [t_id, items] : m_local_items) { + table_ids.push_back(t_id); + } + std::sort(table_ids.begin(), table_ids.end()); + table_ids.erase(std::unique(table_ids.begin(), table_ids.end()), table_ids.end()); + + auto merge_func = [](t_id_vec a, t_id_vec b){ + t_id_vec res; + std::set_union(a.cbegin(), a.cend(), b.cbegin(), b.cend(), std::back_inserter(res)); + return res; + }; + + table_ids = m_comm.all_reduce(table_ids, merge_func); + m_comm.barrier(); + + std::vector local_weights_vec(table_ids.size()); + for (int i = 0; i < table_ids.size(); i++) { + local_weights_vec[i] = m_each_tables_local_weight[table_ids[i]]; + } + MPI_Comm comm = m_comm.get_mpi_comm(); + std::vector w8s_prfx_sum(table_ids.size()); + MPI_Exscan(local_weights_vec.data(), w8s_prfx_sum.data(), local_weights_vec.size(), MPI_DOUBLE, MPI_SUM, comm); + std::vector global_table_w8s(table_ids.size()); // Global total weight for each table + MPI_Allreduce(local_weights_vec.data(), global_table_w8s.data(), local_weights_vec.size(), MPI_DOUBLE, MPI_SUM, comm); + m_comm.barrier(); + // m_comm.cout0("Completed collective communications in balance()"); + + ASSERT_RELEASE(ygm::is_same(table_ids.size(), m_comm)); + ASSERT_RELEASE(ygm::is_same(w8s_prfx_sum.size(), m_comm)); + ASSERT_RELEASE(ygm::is_same(global_table_w8s.size(), m_comm)); + + std::unordered_map> local_table_items; + // Need to initialize all vectors or accessing map inside lambdas fails + for (auto& t_id : table_ids) { + local_table_items.insert({t_id, std::vector()}); + } + using ygm_items_map_ptr = ygm::ygm_ptr>>; + ygm_items_map_ptr ptr_tbl_items = m_comm.make_ygm_ptr(local_table_items); + m_comm.barrier(); // ensure ygm_ptr is created on every rank + + for (uint32_t t = 0; t < table_ids.size(); t++) { + uint32_t t_id = table_ids[t]; + double p = w8s_prfx_sum[t]; + double global_w8 = global_table_w8s[t]; + double target_w8 = global_w8 / m_comm.size(); // target weight per rank + uint32_t dest_rank = p / target_w8; + double curr_weight = std::fmod(p, target_w8); // Spillage weight + std::vector items_to_send; + + std::vector& local_items = m_local_items[t_id]; + for (uint64_t i = 0; i < local_items.size(); i++) { // WARNING: SIZE OF m_local_items CAN GROW DURING LOOP + item local_item = local_items[i]; + if (curr_weight + local_item.weight >= target_w8) { + + double remaining_item_w8 = curr_weight + local_item.weight - target_w8; + double weight_to_send = local_item.weight - remaining_item_w8; + curr_weight += weight_to_send; + item item_to_send = {local_item.t_id, local_item.i_id, weight_to_send}; + items_to_send.push_back(item_to_send); + + // if (dest_rank == m_comm.rank()) { + // local_table_items[t_id].insert(local_table_items[t_id].end(), items_to_send.begin(), items_to_send.end()); + if ((curr_weight > 0.0001) && (dest_rank < m_comm.size())) { // Accounts for rounding errors + // Lambda that moves weights to dest_rank's local_table_items + m_comm.async(dest_rank, [](table_id_type table_id, std::vector items, ygm_items_map_ptr my_items_ptr) { + my_items_ptr->at(table_id).insert(my_items_ptr->at(table_id).end(), items.begin(), items.end()); + }, t_id, items_to_send, ptr_tbl_items); + } + // Handle case where item weight is large enough to span multiple partitions + if (remaining_item_w8 >= target_w8) { + local_items.push_back({local_item.t_id, local_item.i_id, remaining_item_w8}); + curr_weight = 0; + } else { + curr_weight = remaining_item_w8; + } + items_to_send.clear(); + if (curr_weight != 0) { + items_to_send.push_back({local_item.t_id, local_item.i_id, curr_weight}); + } + dest_rank++; + + } else { + items_to_send.push_back(local_item); + curr_weight += local_item.weight; + } + } + + // Need to handle left over items and send them to neccessary destination + if (curr_weight > 0.0001 && dest_rank < m_comm.size()) { // Account for rounding errors + // if (dest_rank == m_comm.rank()) { + // local_table_items[t_id].insert(local_table_items[t_id].end(), items_to_send.begin(), items_to_send.end()); + // } else { + // Lambda that moves the weights to curr_ranks local_table_items + m_comm.async(dest_rank, [](table_id_type table_id, std::vector items, ygm_items_map_ptr my_items_ptr) { + my_items_ptr->at(table_id).insert(my_items_ptr->at(table_id).end(), items.begin(), items.end()); + }, t_id, items_to_send, ptr_tbl_items); + // } + } + } + m_comm.barrier(); + std::swap(local_table_items, m_local_items); + + // Update the total weight of each table that the rank owns + for (uint32_t t = 0; t < table_ids.size(); t++) { + uint32_t t_id = table_ids[t]; + double new_total_w8 = 0; + for (auto& itm : m_local_items[t_id]) { + new_total_w8 += itm.weight; + } + m_each_tables_local_weight[t_id] = new_total_w8; + } + + ASSERT_RELEASE(m_each_tables_local_weight.size() == m_local_items.size()); + m_balanced = true; + } + + bool check_balancing() { + // Not a perfect check. Will check if every rank has the same weight for each table + // but does not confirm that the actual target weight is actually on each rank + std::vector w8_vec; + for (const auto& w : m_each_tables_local_weight) { + w8_vec.push_back(w.second); + // m_comm.cout(w.first, ", ", w.second); + } + + m_comm.barrier(); + + ASSERT_RELEASE(ygm::is_same(w8_vec.size(), m_comm)); + + auto equal_vec = [](std::vector a, std::vector b){ + if (a.size() != b.size()) { + return false; + } else { + for (size_t i = 0; i < a.size(); i++) { + if (std::abs(a[i] - b[i]) > 0.0001) + return false; + } + } + return true; + }; + + bool balanced = ygm::is_same(w8_vec, m_comm, equal_vec); + + return balanced; + } + + void build_alias_tables() { + if (!m_balanced) { + balance_weight(); + } + m_comm.barrier(); + for (auto& [t_id, items] : m_local_items) { + // Make average weight of items 1, so coin flip is simple across all tables + double avg_w8 = m_each_tables_local_weight[t_id] / items.size(); + double new_total_weight = 0; + for (auto& i : items) { + i.weight /= avg_w8; + new_total_weight += i.weight; + } + double new_avg_w8 = new_total_weight / items.size(); // Should be 1 + if (std::abs((new_avg_w8 - 1.0)) < 0.0001) { + ASSERT_RELEASE(std::abs((new_avg_w8 - 1.0)) < 0.0001); + } + + // Implementation of Vose's algorithm, utilized Keith Schwarz numerically stable version + // https://www.keithschwarz.com/darts-dice-coins/ + std::vector heavy_items; + std::vector light_items; + for (auto& item : items) { + if (item.weight < 1.0) { + light_items.push_back(item); + } else { + heavy_items.push_back(item); + } + } + + while (!light_items.empty() && !heavy_items.empty()) { + item l = light_items.back(); + item h = heavy_items.back(); + table_item tbl_i = {l.weight, l.i_id, h.i_id}; + m_local_alias_tables[t_id].push_back(tbl_i); + h.weight = (h.weight + l.weight) - 1; + light_items.pop_back(); + if (h.weight < 1) { + light_items.push_back(h); + heavy_items.pop_back(); + } + } + + // Either heavy items or light_items is empty, need to flush the non empty + // vector and add them to the alias table with a p value of 1 + while (!heavy_items.empty()) { + item h = heavy_items.back(); + table_item tbl_i = {1.0, h.i_id, 0}; + m_local_alias_tables[t_id].push_back(tbl_i); + heavy_items.pop_back(); + } + while (!light_items.empty()) { + item l = light_items.back(); + table_item tbl_i = {1.0, l.i_id, 0}; + m_local_alias_tables[t_id].push_back(tbl_i); + light_items.pop_back(); + } + } + m_comm.barrier(); + m_tables_built = true; + } + + uint32_t get_rank() { + return m_comm.rank(); + } + + template + void async_sample(table_id_type table_id, Visitor visitor, const VisitorArgs &...args) { // Sample should be provided a function with takes as an argument the item type + + ASSERT_RELEASE(m_tables_built); + + auto visit_wrapper = [](auto ptr_MAT, table_id_type t_id, const VisitorArgs &...args) { + std::uniform_int_distribution table_item_dist(0, ptr_MAT->m_local_alias_tables[t_id].size()-1); + table_item tbl_itm = ptr_MAT->m_local_alias_tables[t_id][table_item_dist(ptr_MAT->m_rng)]; + item_id_type s; + if (tbl_itm.p == 1) { + s = tbl_itm.a; + } else { + std::uniform_real_distribution zero_one_dist(0.0, 1.0); + float f = zero_one_dist(ptr_MAT->m_rng); + if (f < tbl_itm.p) { + s = tbl_itm.a; + } else { + s = tbl_itm.b; + } + } + Visitor *vis = nullptr; + ygm::meta::apply_optional(*vis, std::make_tuple(ptr_MAT), std::forward_as_tuple(s, t_id, args...)); + }; + + uint32_t dest_rank = m_rank_dist(m_rng); + m_comm.async(dest_rank, visit_wrapper, pthis, table_id, std::forward(args)...); + } + + private: + bool m_balanced; + bool m_tables_built; + ygm::comm& m_comm; + ygm::ygm_ptr pthis; + std::uniform_int_distribution m_rank_dist; + RNG m_rng; + std::unordered_map> m_local_items; + std::unordered_map> m_local_alias_tables; + std::map m_each_tables_local_weight; +}; + +} // namespace ygm::random \ No newline at end of file diff --git a/include/ygm/random.hpp b/include/ygm/random/random.hpp similarity index 95% rename from include/ygm/random.hpp rename to include/ygm/random/random.hpp index 8edaeacd..5efe13c0 100644 --- a/include/ygm/random.hpp +++ b/include/ygm/random/random.hpp @@ -7,7 +7,7 @@ #include -namespace ygm { +namespace ygm::random { /// @brief A simple offset rng alias /// @tparam RandomEngine The underlying random engine, e.g. std::mt19937 diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt index 0ba0414e..a117575d 100644 --- a/test/CMakeLists.txt +++ b/test/CMakeLists.txt @@ -66,6 +66,7 @@ add_ygm_test(test_recursion_progress) add_ygm_test(test_gather_topk) add_ygm_test(test_reduce) add_ygm_test(test_transform) +add_ygm_test(test_alias_table) if (Boost_FOUND) add_ygm_seq_test(test_cereal_boost_json) diff --git a/test/test_alias_table.cpp b/test/test_alias_table.cpp new file mode 100644 index 00000000..743d30de --- /dev/null +++ b/test/test_alias_table.cpp @@ -0,0 +1,80 @@ +// Copyright 2019-2021 Lawrence Livermore National Security, LLC and other YGM +// Project Developers. See the top-level COPYRIGHT file for details. +// +// SPDX-License-Identifier: MIT + +#undef NDEBUG + +#include +#include +#include +#include +#include + +int main(int argc, char** argv) { + + ygm::comm world(&argc, &argv); + using YGM_RNG = ygm::random::default_random_engine<>; + int seed = 150; + YGM_RNG ygm_rng(world, seed); + { + ygm::container::bag> bag_of_items(world); + + uint32_t n_items_per_rank = 1000; + const int max_item_weight = 100; + + std::uniform_real_distribution dist(0, max_item_weight); + const uint32_t rank = world.rank(); + for (uint32_t i = 0; i < n_items_per_rank; i++) { + uint32_t id = i + rank * n_items_per_rank; + double w = dist(ygm_rng); + bag_of_items.async_insert({id, w}); + } + world.barrier(); + + ygm::random::alias_table alias_tbl(world, ygm_rng, bag_of_items); + + static uint32_t samples; + uint32_t samples_per_rank = 1000; + for (int i = 0; i < samples_per_rank; i++) { + alias_tbl.async_sample([](auto ptr, uint32_t item){ + samples++; + }); + } + world.barrier(); + uint32_t total_samples = ygm::sum(samples, world); + YGM_ASSERT_RELEASE(total_samples == (samples_per_rank * world.size())); + } + + { + ygm::container::map map_of_items(world); + + uint32_t n_items_per_rank = 1000; + const int max_item_weight = 100; + + std::uniform_real_distribution dist(0, max_item_weight); + const uint32_t rank = world.rank(); + for (uint32_t i = 0; i < n_items_per_rank; i++) { + uint32_t id = i + rank * n_items_per_rank; + double w = dist(ygm_rng); + map_of_items.async_insert(id, w); + } + world.barrier(); + + ygm::random::alias_table alias_tbl(world, ygm_rng, map_of_items); + + static uint32_t samples; + uint32_t samples_per_rank = 1000; + for (int i = 0; i < samples_per_rank; i++) { + alias_tbl.async_sample([](auto ptr, uint32_t item){ + samples++; + }); + } + world.barrier(); + uint32_t total_samples = ygm::sum(samples, world); + YGM_ASSERT_RELEASE(total_samples == (samples_per_rank * world.size())); + + } + + return 0; +} diff --git a/test/test_bag.cpp b/test/test_bag.cpp index 04362208..6790c706 100644 --- a/test/test_bag.cpp +++ b/test/test_bag.cpp @@ -10,7 +10,7 @@ #include #include #include -#include +#include int main(int argc, char** argv) { ygm::comm world(&argc, &argv); @@ -170,12 +170,12 @@ int main(int argc, char** argv) { } } int seed = 100; - ygm::default_random_engine<> rng1 = - ygm::default_random_engine<>(world, seed); + ygm::random::default_random_engine<> rng1 = + ygm::random::default_random_engine<>(world, seed); bbag.local_shuffle(rng1); - ygm::default_random_engine<> rng2 = - ygm::default_random_engine<>(world, seed); + ygm::random::default_random_engine<> rng2 = + ygm::random::default_random_engine<>(world, seed); bbag.global_shuffle(rng2); bbag.local_shuffle(); diff --git a/test/test_gather_topk.cpp b/test/test_gather_topk.cpp index 99a5a4dd..a65caff0 100644 --- a/test/test_gather_topk.cpp +++ b/test/test_gather_topk.cpp @@ -11,7 +11,6 @@ #include #include #include -#include int main(int argc, char** argv) { ygm::comm world(&argc, &argv); diff --git a/test/test_random.cpp b/test/test_random.cpp index 9f28e714..cf926d10 100644 --- a/test/test_random.cpp +++ b/test/test_random.cpp @@ -6,7 +6,7 @@ #undef NDEBUG #include #include -#include +#include #include @@ -17,7 +17,7 @@ int main(int argc, char** argv) { // Test default_random_engine { std::uint32_t seed = 100; - ygm::default_random_engine<> rng(world, seed); + ygm::random::default_random_engine<> rng(world, seed); ygm::container::counting_set seed_set(world); ygm::container::counting_set rn_set(world); ygm::container::counting_set sample_set(world); diff --git a/test/test_reduce.cpp b/test/test_reduce.cpp index c43861cb..9ca06645 100644 --- a/test/test_reduce.cpp +++ b/test/test_reduce.cpp @@ -11,7 +11,6 @@ #include #include #include -#include int main(int argc, char** argv) { ygm::comm world(&argc, &argv); diff --git a/test/test_transform.cpp b/test/test_transform.cpp index 607916d9..2fd6d6c6 100644 --- a/test/test_transform.cpp +++ b/test/test_transform.cpp @@ -11,7 +11,6 @@ #include #include #include -#include int main(int argc, char** argv) { ygm::comm world(&argc, &argv); From e6612b2a73fbe7b0cebb59fd8c303e3b0c0d105d Mon Sep 17 00:00:00 2001 From: Lance Fletcher Date: Wed, 2 Jul 2025 22:41:01 -0700 Subject: [PATCH 02/15] Modified types names and cleaned up a few other things. Still need to clean up the use of concepts. --- .gitignore | 2 + include/ygm/random/alias_table.hpp | 88 +++--- include/ygm/random/multi_alias_table.hpp | 364 ----------------------- test/test_alias_table.cpp | 7 +- 4 files changed, 51 insertions(+), 410 deletions(-) delete mode 100644 include/ygm/random/multi_alias_table.hpp diff --git a/.gitignore b/.gitignore index cb84f20f..15752655 100644 --- a/.gitignore +++ b/.gitignore @@ -4,3 +4,5 @@ build* .vscode* .idea* .cache/ +CMakeCache.txt +CMakeFiles/ diff --git a/include/ygm/random/alias_table.hpp b/include/ygm/random/alias_table.hpp index 4818ccf3..ca873f95 100644 --- a/include/ygm/random/alias_table.hpp +++ b/include/ygm/random/alias_table.hpp @@ -1,7 +1,7 @@ #include #include -#include -#include +#include +#include #include #include #include @@ -11,17 +11,9 @@ namespace ygm::random { - - template concept convertable_to_double = std::convertible_to; -template -concept is_pair = requires { - typename T::first_type; - typename T::second_type; -} && std::same_as>; - template concept is_associative_container = requires { typename T::key_type; @@ -33,18 +25,34 @@ concept is_value_container = requires { typename T::value_type; }; +// template +// concept valid_stl_associative_container = requires { +// typename T::key_type; +// typename T::value_type; +// }; + + +// template +// concept valid_stl_value_container = requires { +// typename T::value_type; +// requires ygm::container::detail::DoubleItemTuple; +// }; + +// template +// concept valid_stl_container = +// is_associative_container || valid_stl_value_container; -template +template class alias_table { - using self_type = alias_table; - using weight_type = double; + using self_type = alias_table; + using weight_t = double; public: struct item { - item_type id; - weight_type weight; + Item id; + weight_t weight; template void serialize(Archive& ar) { // Needed for cereal serialization @@ -54,38 +62,35 @@ class alias_table { struct table_item { double p; // prob item a is selected = p, prob item b is selected = (1 - p) - item_type a; - item_type b; - - template - void serialize(Archive& ar) { // Needed for cereal serialization - ar(p, a, b); - } + Item a; + Item b; }; template - requires is_value_container && - is_pair && - convertable_to_double + requires ygm::container::detail::HasForAll && + is_value_container && + ygm::container::detail::DoubleItemTuple && + convertable_to_double> alias_table(ygm::comm &comm, RNG &rng, Container &c) : m_comm(comm), pthis(this), m_rng(rng), m_rank_dist(0, comm.size()-1) { c.for_all([&](const auto& id_weight_pair){ - m_local_items.push_back({id_weight_pair.first, id_weight_pair.second}); + m_local_items.push_back({std::get<0>(id_weight_pair), std::get<1>(id_weight_pair)}); }); comm.barrier(); balance_weight(); comm.barrier(); - YGM_ASSERT_RELEASE(check_balancing()); + YGM_ASSERT_RELEASE(is_balanced()); build_alias_table(); m_local_items.clear(); } template - requires is_associative_container && + requires ygm::container::detail::HasForAll && + is_associative_container && convertable_to_double alias_table(ygm::comm &comm, RNG &rng, Container &c) : m_comm(comm), @@ -98,19 +103,12 @@ class alias_table { comm.barrier(); balance_weight(); comm.barrier(); - YGM_ASSERT_RELEASE(check_balancing()); + YGM_ASSERT_RELEASE(is_balanced()); build_alias_table(); m_local_items.clear(); } - // For debugging - std::string get_weight_str(std::vector items) { - std::string weights_str = "< "; - for (auto& itm : items) { - weights_str += std::to_string(itm.weight) + " "; - } - return weights_str + ">"; - } + private: void balance_weight() { @@ -183,7 +181,7 @@ class alias_table { std::swap(new_local_items, m_local_items); } - bool check_balancing() { + bool is_balanced() { // Not a perfect check. Will check if every rank has the same weight for each table // but does not confirm that the actual target weight is met on each rank double local_weight = 0.0; @@ -257,19 +255,21 @@ class alias_table { } m_comm.barrier(); } + + public: template void async_sample(Visitor visitor, const VisitorArgs &...args) { // Sample should be provided a lambda/functor that takes as an argument the item type - auto sample_wrapper = [](auto ptr_MAT, const VisitorArgs &...args) { - std::uniform_int_distribution table_item_dist(0, ptr_MAT->m_local_alias_table.size()-1); - table_item tbl_itm = ptr_MAT->m_local_alias_table[table_item_dist(ptr_MAT->m_rng)]; - item_type s; + auto sample_wrapper = [](auto ptr_a_tbl, const VisitorArgs &...args) { + std::uniform_int_distribution table_item_dist(0, ptr_a_tbl->m_local_alias_table.size()-1); + table_item tbl_itm = ptr_a_tbl->m_local_alias_table[table_item_dist(ptr_a_tbl->m_rng)]; + Item s; if (tbl_itm.p == 1) { s = tbl_itm.a; } else { std::uniform_real_distribution zero_one_dist(0.0, 1.0); - float f = zero_one_dist(ptr_MAT->m_rng); + float f = zero_one_dist(ptr_a_tbl->m_rng); if (f < tbl_itm.p) { s = tbl_itm.a; } else { @@ -277,7 +277,7 @@ class alias_table { } } Visitor *vis = nullptr; - ygm::meta::apply_optional(*vis, std::make_tuple(ptr_MAT), std::forward_as_tuple(s, args...)); + ygm::meta::apply_optional(*vis, std::make_tuple(ptr_a_tbl), std::forward_as_tuple(s, args...)); }; uint32_t dest_rank = m_rank_dist(m_rng); diff --git a/include/ygm/random/multi_alias_table.hpp b/include/ygm/random/multi_alias_table.hpp deleted file mode 100644 index 4faeb92d..00000000 --- a/include/ygm/random/multi_alias_table.hpp +++ /dev/null @@ -1,364 +0,0 @@ -#include -#include -#include -#include -#include -#include -#include -#include -#include - - -namespace ygm::random { - -template -class multi_alias_table { - - using self_type = multi_alias_table; - - public: - - struct item { - table_id_type t_id; - item_id_type i_id; - weight_type weight; - - template - void serialize(Archive& ar) { // Needed for cereal serialization - ar(t_id, i_id, weight); - } - }; - - - struct table_item { - double p; // prob item a is selected = p, prob item b is selected = (1 - p) - item_id_type a; - item_id_type b; - - template - void serialize(Archive& ar) { // Needed for cereal serialization - ar(p, a, b); - } - }; - - - multi_alias_table(ygm::comm &comm, RNG rng) : - m_comm(comm), - pthis(this), - m_rng(rng), - m_tables_built(false), - m_balanced(false), - // m_local_weight(0), - m_rank_dist(0, comm.size()-1) { } - - - template - void load_items_from_container(container &cont) { - cont.for_all([&](item itm){ - m_local_items[itm.t_id].push_back(itm); - m_each_tables_local_weight[itm.t_id] += itm.weight; - }); - } - - void local_add_item (std::tuple i) { - item itm({std::get<0>(i), std::get<1>(i), std::get<2>(i)}); - m_local_items[itm.t_id].push_back(itm); - m_each_tables_local_weight[itm.t_id] += itm.weight; - } - - void add_item(std::tuple i) { - item itm({std::get<0>(i), std::get<1>(i), std::get<2>(i)}); - - auto add_item_to_rand_rank = [](item i, auto ptr_MAT){ - ptr_MAT->m_local_items[i.t_id].push_back(i); - ptr_MAT->m_each_tables_local_weight[i.t_id] += i.weight; - }; - - m_comm.async(m_rank_dist(m_rng), add_item_to_rand_rank, itm, pthis); - } - - // For debugging - std::string get_weight_str(std::vector items) { - std::string weights_str = "< "; - for (auto& itm : items) { - weights_str += std::to_string(itm.weight) + " "; - } - return weights_str + ">"; - } - - std::vector get_table_ids() { - std::vector table_ids; - for (const table_id_type t : m_each_tables_local_weight) { - table_ids.push_back(t); - } - return table_ids; - } - - void clear() { - m_comm.barrier(); - m_tables_built = false; - m_local_items.clear(); - m_local_alias_tables.clear(); - m_each_tables_local_weight.clear(); - } - - void clear_items() { - m_local_items.clear(); - } - - void balance_weight() { - using t_id_vec = std::vector; - t_id_vec table_ids; - for (const auto& [t_id, items] : m_local_items) { - table_ids.push_back(t_id); - } - std::sort(table_ids.begin(), table_ids.end()); - table_ids.erase(std::unique(table_ids.begin(), table_ids.end()), table_ids.end()); - - auto merge_func = [](t_id_vec a, t_id_vec b){ - t_id_vec res; - std::set_union(a.cbegin(), a.cend(), b.cbegin(), b.cend(), std::back_inserter(res)); - return res; - }; - - table_ids = m_comm.all_reduce(table_ids, merge_func); - m_comm.barrier(); - - std::vector local_weights_vec(table_ids.size()); - for (int i = 0; i < table_ids.size(); i++) { - local_weights_vec[i] = m_each_tables_local_weight[table_ids[i]]; - } - MPI_Comm comm = m_comm.get_mpi_comm(); - std::vector w8s_prfx_sum(table_ids.size()); - MPI_Exscan(local_weights_vec.data(), w8s_prfx_sum.data(), local_weights_vec.size(), MPI_DOUBLE, MPI_SUM, comm); - std::vector global_table_w8s(table_ids.size()); // Global total weight for each table - MPI_Allreduce(local_weights_vec.data(), global_table_w8s.data(), local_weights_vec.size(), MPI_DOUBLE, MPI_SUM, comm); - m_comm.barrier(); - // m_comm.cout0("Completed collective communications in balance()"); - - ASSERT_RELEASE(ygm::is_same(table_ids.size(), m_comm)); - ASSERT_RELEASE(ygm::is_same(w8s_prfx_sum.size(), m_comm)); - ASSERT_RELEASE(ygm::is_same(global_table_w8s.size(), m_comm)); - - std::unordered_map> local_table_items; - // Need to initialize all vectors or accessing map inside lambdas fails - for (auto& t_id : table_ids) { - local_table_items.insert({t_id, std::vector()}); - } - using ygm_items_map_ptr = ygm::ygm_ptr>>; - ygm_items_map_ptr ptr_tbl_items = m_comm.make_ygm_ptr(local_table_items); - m_comm.barrier(); // ensure ygm_ptr is created on every rank - - for (uint32_t t = 0; t < table_ids.size(); t++) { - uint32_t t_id = table_ids[t]; - double p = w8s_prfx_sum[t]; - double global_w8 = global_table_w8s[t]; - double target_w8 = global_w8 / m_comm.size(); // target weight per rank - uint32_t dest_rank = p / target_w8; - double curr_weight = std::fmod(p, target_w8); // Spillage weight - std::vector items_to_send; - - std::vector& local_items = m_local_items[t_id]; - for (uint64_t i = 0; i < local_items.size(); i++) { // WARNING: SIZE OF m_local_items CAN GROW DURING LOOP - item local_item = local_items[i]; - if (curr_weight + local_item.weight >= target_w8) { - - double remaining_item_w8 = curr_weight + local_item.weight - target_w8; - double weight_to_send = local_item.weight - remaining_item_w8; - curr_weight += weight_to_send; - item item_to_send = {local_item.t_id, local_item.i_id, weight_to_send}; - items_to_send.push_back(item_to_send); - - // if (dest_rank == m_comm.rank()) { - // local_table_items[t_id].insert(local_table_items[t_id].end(), items_to_send.begin(), items_to_send.end()); - if ((curr_weight > 0.0001) && (dest_rank < m_comm.size())) { // Accounts for rounding errors - // Lambda that moves weights to dest_rank's local_table_items - m_comm.async(dest_rank, [](table_id_type table_id, std::vector items, ygm_items_map_ptr my_items_ptr) { - my_items_ptr->at(table_id).insert(my_items_ptr->at(table_id).end(), items.begin(), items.end()); - }, t_id, items_to_send, ptr_tbl_items); - } - // Handle case where item weight is large enough to span multiple partitions - if (remaining_item_w8 >= target_w8) { - local_items.push_back({local_item.t_id, local_item.i_id, remaining_item_w8}); - curr_weight = 0; - } else { - curr_weight = remaining_item_w8; - } - items_to_send.clear(); - if (curr_weight != 0) { - items_to_send.push_back({local_item.t_id, local_item.i_id, curr_weight}); - } - dest_rank++; - - } else { - items_to_send.push_back(local_item); - curr_weight += local_item.weight; - } - } - - // Need to handle left over items and send them to neccessary destination - if (curr_weight > 0.0001 && dest_rank < m_comm.size()) { // Account for rounding errors - // if (dest_rank == m_comm.rank()) { - // local_table_items[t_id].insert(local_table_items[t_id].end(), items_to_send.begin(), items_to_send.end()); - // } else { - // Lambda that moves the weights to curr_ranks local_table_items - m_comm.async(dest_rank, [](table_id_type table_id, std::vector items, ygm_items_map_ptr my_items_ptr) { - my_items_ptr->at(table_id).insert(my_items_ptr->at(table_id).end(), items.begin(), items.end()); - }, t_id, items_to_send, ptr_tbl_items); - // } - } - } - m_comm.barrier(); - std::swap(local_table_items, m_local_items); - - // Update the total weight of each table that the rank owns - for (uint32_t t = 0; t < table_ids.size(); t++) { - uint32_t t_id = table_ids[t]; - double new_total_w8 = 0; - for (auto& itm : m_local_items[t_id]) { - new_total_w8 += itm.weight; - } - m_each_tables_local_weight[t_id] = new_total_w8; - } - - ASSERT_RELEASE(m_each_tables_local_weight.size() == m_local_items.size()); - m_balanced = true; - } - - bool check_balancing() { - // Not a perfect check. Will check if every rank has the same weight for each table - // but does not confirm that the actual target weight is actually on each rank - std::vector w8_vec; - for (const auto& w : m_each_tables_local_weight) { - w8_vec.push_back(w.second); - // m_comm.cout(w.first, ", ", w.second); - } - - m_comm.barrier(); - - ASSERT_RELEASE(ygm::is_same(w8_vec.size(), m_comm)); - - auto equal_vec = [](std::vector a, std::vector b){ - if (a.size() != b.size()) { - return false; - } else { - for (size_t i = 0; i < a.size(); i++) { - if (std::abs(a[i] - b[i]) > 0.0001) - return false; - } - } - return true; - }; - - bool balanced = ygm::is_same(w8_vec, m_comm, equal_vec); - - return balanced; - } - - void build_alias_tables() { - if (!m_balanced) { - balance_weight(); - } - m_comm.barrier(); - for (auto& [t_id, items] : m_local_items) { - // Make average weight of items 1, so coin flip is simple across all tables - double avg_w8 = m_each_tables_local_weight[t_id] / items.size(); - double new_total_weight = 0; - for (auto& i : items) { - i.weight /= avg_w8; - new_total_weight += i.weight; - } - double new_avg_w8 = new_total_weight / items.size(); // Should be 1 - if (std::abs((new_avg_w8 - 1.0)) < 0.0001) { - ASSERT_RELEASE(std::abs((new_avg_w8 - 1.0)) < 0.0001); - } - - // Implementation of Vose's algorithm, utilized Keith Schwarz numerically stable version - // https://www.keithschwarz.com/darts-dice-coins/ - std::vector heavy_items; - std::vector light_items; - for (auto& item : items) { - if (item.weight < 1.0) { - light_items.push_back(item); - } else { - heavy_items.push_back(item); - } - } - - while (!light_items.empty() && !heavy_items.empty()) { - item l = light_items.back(); - item h = heavy_items.back(); - table_item tbl_i = {l.weight, l.i_id, h.i_id}; - m_local_alias_tables[t_id].push_back(tbl_i); - h.weight = (h.weight + l.weight) - 1; - light_items.pop_back(); - if (h.weight < 1) { - light_items.push_back(h); - heavy_items.pop_back(); - } - } - - // Either heavy items or light_items is empty, need to flush the non empty - // vector and add them to the alias table with a p value of 1 - while (!heavy_items.empty()) { - item h = heavy_items.back(); - table_item tbl_i = {1.0, h.i_id, 0}; - m_local_alias_tables[t_id].push_back(tbl_i); - heavy_items.pop_back(); - } - while (!light_items.empty()) { - item l = light_items.back(); - table_item tbl_i = {1.0, l.i_id, 0}; - m_local_alias_tables[t_id].push_back(tbl_i); - light_items.pop_back(); - } - } - m_comm.barrier(); - m_tables_built = true; - } - - uint32_t get_rank() { - return m_comm.rank(); - } - - template - void async_sample(table_id_type table_id, Visitor visitor, const VisitorArgs &...args) { // Sample should be provided a function with takes as an argument the item type - - ASSERT_RELEASE(m_tables_built); - - auto visit_wrapper = [](auto ptr_MAT, table_id_type t_id, const VisitorArgs &...args) { - std::uniform_int_distribution table_item_dist(0, ptr_MAT->m_local_alias_tables[t_id].size()-1); - table_item tbl_itm = ptr_MAT->m_local_alias_tables[t_id][table_item_dist(ptr_MAT->m_rng)]; - item_id_type s; - if (tbl_itm.p == 1) { - s = tbl_itm.a; - } else { - std::uniform_real_distribution zero_one_dist(0.0, 1.0); - float f = zero_one_dist(ptr_MAT->m_rng); - if (f < tbl_itm.p) { - s = tbl_itm.a; - } else { - s = tbl_itm.b; - } - } - Visitor *vis = nullptr; - ygm::meta::apply_optional(*vis, std::make_tuple(ptr_MAT), std::forward_as_tuple(s, t_id, args...)); - }; - - uint32_t dest_rank = m_rank_dist(m_rng); - m_comm.async(dest_rank, visit_wrapper, pthis, table_id, std::forward(args)...); - } - - private: - bool m_balanced; - bool m_tables_built; - ygm::comm& m_comm; - ygm::ygm_ptr pthis; - std::uniform_int_distribution m_rank_dist; - RNG m_rng; - std::unordered_map> m_local_items; - std::unordered_map> m_local_alias_tables; - std::map m_each_tables_local_weight; -}; - -} // namespace ygm::random \ No newline at end of file diff --git a/test/test_alias_table.cpp b/test/test_alias_table.cpp index 743d30de..9617aff7 100644 --- a/test/test_alias_table.cpp +++ b/test/test_alias_table.cpp @@ -15,10 +15,13 @@ int main(int argc, char** argv) { ygm::comm world(&argc, &argv); using YGM_RNG = ygm::random::default_random_engine<>; - int seed = 150; + int seed = 42; YGM_RNG ygm_rng(world, seed); + + // Test constructor works with YGM value container { ygm::container::bag> bag_of_items(world); + // ygm::container::bag> bag_of_items(world); uint32_t n_items_per_rank = 1000; const int max_item_weight = 100; @@ -46,6 +49,7 @@ int main(int argc, char** argv) { YGM_ASSERT_RELEASE(total_samples == (samples_per_rank * world.size())); } + // Test constructor works with YGM associative container { ygm::container::map map_of_items(world); @@ -73,7 +77,6 @@ int main(int argc, char** argv) { world.barrier(); uint32_t total_samples = ygm::sum(samples, world); YGM_ASSERT_RELEASE(total_samples == (samples_per_rank * world.size())); - } return 0; From 8f55fca67843e7b8bb4f591d83f966e7c13e1501 Mon Sep 17 00:00:00 2001 From: Lance Fletcher Date: Fri, 11 Jul 2025 21:42:08 -0700 Subject: [PATCH 03/15] Modified files relying on ygm/random.hpp to point to ygm/random/random.hpp --- test/test_transform.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/test/test_transform.cpp b/test/test_transform.cpp index a8979c13..ac70371e 100644 --- a/test/test_transform.cpp +++ b/test/test_transform.cpp @@ -12,7 +12,7 @@ #include #include #include -#include +#include int main(int argc, char** argv) { ygm::comm world(&argc, &argv); From 4f9f22a76279c859e70061f3a4e1d535eec71e31 Mon Sep 17 00:00:00 2001 From: Lance Fletcher Date: Wed, 10 Dec 2025 14:39:45 -0800 Subject: [PATCH 04/15] Changed 4 space indentation to 2 space. Also moved ygm::detail::random_engine to ygm::random::random_engine --- include/ygm/random/alias_table.hpp | 485 +++++++++++++++-------------- include/ygm/random/random.hpp | 43 ++- test/test_alias_table.cpp | 13 +- 3 files changed, 293 insertions(+), 248 deletions(-) diff --git a/include/ygm/random/alias_table.hpp b/include/ygm/random/alias_table.hpp index ca873f95..2113fbf1 100644 --- a/include/ygm/random/alias_table.hpp +++ b/include/ygm/random/alias_table.hpp @@ -1,3 +1,11 @@ + +// Copyright 2019-2025 Lawrence Livermore National Security, LLC and other YGM +// Project Developers. See the top-level COPYRIGHT file for details. +// +// SPDX-License-Identifier: MIT + +#pragma once + #include #include #include @@ -45,252 +53,247 @@ concept is_value_container = requires { template class alias_table { - using self_type = alias_table; - using weight_t = double; - - public: - - struct item { - Item id; - weight_t weight; - - template - void serialize(Archive& ar) { // Needed for cereal serialization - ar(id, weight); - } - }; - - struct table_item { - double p; // prob item a is selected = p, prob item b is selected = (1 - p) - Item a; - Item b; - }; - - - template - requires ygm::container::detail::HasForAll && - is_value_container && - ygm::container::detail::DoubleItemTuple && - convertable_to_double> - alias_table(ygm::comm &comm, RNG &rng, Container &c) : - m_comm(comm), - pthis(this), - m_rng(rng), - m_rank_dist(0, comm.size()-1) { - c.for_all([&](const auto& id_weight_pair){ - m_local_items.push_back({std::get<0>(id_weight_pair), std::get<1>(id_weight_pair)}); - }); - comm.barrier(); - balance_weight(); - comm.barrier(); - YGM_ASSERT_RELEASE(is_balanced()); - build_alias_table(); - m_local_items.clear(); + public: + using self_type = alias_table; + using weight_t = double; + + struct item { + Item id; + weight_t weight; + + template + void serialize(Archive& ar) { // Needed for cereal serialization + ar(id, weight); + } + }; + + struct table_item { + double p; // prob item a is selected = p, prob item b is selected = (1 - p) + Item a; + Item b; + }; + + template + requires ygm::container::detail::HasForAll && + is_value_container && + ygm::container::detail::DoubleItemTuple && + convertable_to_double> + alias_table(ygm::comm &comm, RNG &rng, Container &c) + : m_comm(comm), pthis(this), m_rng(rng), + m_rank_dist(0, comm.size()-1), m_zero_one_dist(0.0, 1.0) { + + c.for_all([&](const auto& id_weight_pair){ + m_local_items.push_back({std::get<0>(id_weight_pair), std::get<1>(id_weight_pair)}); + }); + comm.barrier(); + balance_weight(); + comm.barrier(); + YGM_ASSERT_RELEASE(is_balanced()); + build_alias_table(); + m_local_items.clear(); + } + + template + requires ygm::container::detail::HasForAll && + is_associative_container && + convertable_to_double + alias_table(ygm::comm &comm, RNG &rng, Container &c) + : m_comm(comm), pthis(this), m_rng(rng), + m_rank_dist(0, comm.size()-1), m_zero_one_dist(0.0, 1.0) { + + c.for_all([&](const auto &key, const auto &value){ + m_local_items.push_back({key, value}); + }); + comm.barrier(); + balance_weight(); + comm.barrier(); + YGM_ASSERT_RELEASE(is_balanced()); + build_alias_table(); + m_local_items.clear(); + } + + private: + + void balance_weight() { + MPI_Comm comm = m_comm.get_mpi_comm(); + double local_weight = 0.0; + for (uint32_t i = 0; i < m_local_items.size(); i++) { + local_weight += m_local_items[i].weight; + } + + double global_weight; + MPI_Allreduce(&local_weight, &global_weight, 1, MPI_DOUBLE, MPI_SUM, comm); + + double prfx_sum_weight; + MPI_Exscan(&local_weight, &prfx_sum_weight, 1, MPI_DOUBLE, MPI_SUM, comm); + + + double target_weight = global_weight / m_comm.size(); // target weight per rank + uint32_t dest_rank = prfx_sum_weight / target_weight; + double curr_weight = std::fmod(prfx_sum_weight, target_weight); // Spillage weight + + std::vector new_local_items; + using ygm_items_ptr = ygm::ygm_ptr>; + ygm_items_ptr ptr_new_items = m_comm.make_ygm_ptr(new_local_items); + m_comm.barrier(); // ensure ygm_ptr is created on every rank + + std::vector items_to_send; + for (uint64_t i = 0; i < m_local_items.size(); i++) { // WARNING: size of m_local_items can grow during loop + item local_item = m_local_items[i]; + if (curr_weight + local_item.weight >= target_weight) { + double remaining_weight = curr_weight + local_item.weight - target_weight; + double weight_to_send = local_item.weight - remaining_weight; + curr_weight += weight_to_send; + item item_to_send = {local_item.id, weight_to_send}; + items_to_send.push_back(item_to_send); + + if ((curr_weight > 0.0001) && (dest_rank < m_comm.size())) { // Accounts for rounding errors + // Moves weights to dest_rank's new_local_items + m_comm.async(dest_rank, [](std::vector items, ygm_items_ptr new_items_ptr) { + new_items_ptr->insert(new_items_ptr->end(), items.begin(), items.end()); + }, items_to_send, ptr_new_items); } - template - requires ygm::container::detail::HasForAll && - is_associative_container && - convertable_to_double - alias_table(ygm::comm &comm, RNG &rng, Container &c) : - m_comm(comm), - pthis(this), - m_rng(rng), - m_rank_dist(0, comm.size()-1) { - c.for_all([&](const auto &key, const auto &value){ - m_local_items.push_back({key, value}); - }); - comm.barrier(); - balance_weight(); - comm.barrier(); - YGM_ASSERT_RELEASE(is_balanced()); - build_alias_table(); - m_local_items.clear(); + // Handle case where item weight is large enough to span multiple rank's alias tables + if (remaining_weight >= target_weight) { + m_local_items.push_back({local_item.id, remaining_weight}); + curr_weight = 0; + } else { + curr_weight = remaining_weight; } - - private: - - void balance_weight() { - - MPI_Comm comm = m_comm.get_mpi_comm(); - double local_weight = 0.0; - for (int i = 0; i < m_local_items.size(); i++) { - local_weight += m_local_items[i].weight; - } - - double global_weight; - MPI_Allreduce(&local_weight, &global_weight, 1, MPI_DOUBLE, MPI_SUM, comm); - - double prfx_sum_weight; - MPI_Exscan(&local_weight, &prfx_sum_weight, 1, MPI_DOUBLE, MPI_SUM, comm); - - - double target_weight = global_weight / m_comm.size(); // target weight per rank - uint32_t dest_rank = prfx_sum_weight / target_weight; - double curr_weight = std::fmod(prfx_sum_weight, target_weight); // Spillage weight - - std::vector new_local_items; - using ygm_items_ptr = ygm::ygm_ptr>; - ygm_items_ptr ptr_new_items = m_comm.make_ygm_ptr(new_local_items); - m_comm.barrier(); // ensure ygm_ptr is created on every rank - - std::vector items_to_send; - for (uint64_t i = 0; i < m_local_items.size(); i++) { // WARNING: size of m_local_items can grow during loop - item local_item = m_local_items[i]; - if (curr_weight + local_item.weight >= target_weight) { - - double remaining_weight = curr_weight + local_item.weight - target_weight; - double weight_to_send = local_item.weight - remaining_weight; - curr_weight += weight_to_send; - item item_to_send = {local_item.id, weight_to_send}; - items_to_send.push_back(item_to_send); - - if ((curr_weight > 0.0001) && (dest_rank < m_comm.size())) { // Accounts for rounding errors - // Moves weights to dest_rank's new_local_items - m_comm.async(dest_rank, [](std::vector items, ygm_items_ptr new_items_ptr) { - new_items_ptr->insert(new_items_ptr->end(), items.begin(), items.end()); - }, items_to_send, ptr_new_items); - } - // Handle case where item weight is large enough to span multiple rank's alias tables - if (remaining_weight >= target_weight) { - m_local_items.push_back({local_item.id, remaining_weight}); - curr_weight = 0; - } else { - curr_weight = remaining_weight; - } - items_to_send.clear(); - if (curr_weight != 0) { - items_to_send.push_back({local_item.id, curr_weight}); - } - dest_rank++; - - } else { - items_to_send.push_back(local_item); - curr_weight += local_item.weight; - } - } - - // Need to handle items left in items to send - if (items_to_send.size() > 0 && curr_weight > 0.00001 && dest_rank < m_comm.size()) { // Account for floating point errors - m_comm.async(dest_rank, [](std::vector items, ygm_items_ptr new_items_ptr) { - new_items_ptr->insert(new_items_ptr->end(), items.begin(), items.end()); - }, items_to_send, ptr_new_items); - } - - m_comm.barrier(); - std::swap(new_local_items, m_local_items); - } - - bool is_balanced() { - // Not a perfect check. Will check if every rank has the same weight for each table - // but does not confirm that the actual target weight is met on each rank - double local_weight = 0.0; - for (const auto& itm : m_local_items) { - local_weight += itm.weight; - } - m_comm.barrier(); - auto equal = [](double a, double b){ - return (std::abs(a - b) < 0.00001); - }; - bool balanced = ygm::is_same(local_weight, m_comm, equal); - return balanced; - } - - void build_alias_table() { - // Make average weight of items 1, so coin flip is simple across all tables - double local_weight = 0.0; - for (const auto& itm : m_local_items) { - local_weight += itm.weight; - } - double avg_weight = local_weight / m_local_items.size(); - - double new_total_weight = 0; - for (auto& itm : m_local_items) { - itm.weight /= avg_weight; - new_total_weight += itm.weight; - } - double new_avg_weight = new_total_weight / m_local_items.size(); // Should be 1 - if (std::abs((new_avg_weight - 1.0)) < 0.0001) { - YGM_ASSERT_RELEASE(std::abs((new_avg_weight - 1.0)) < 0.00001); - } - - // Implementation of Vose's algorithm, utilized Keith Schwarz numerically stable version - // https://www.keithschwarz.com/darts-dice-coins/ - std::vector heavy_items; - std::vector light_items; - for (auto& itm : m_local_items) { - if (itm.weight < 1.0) { - light_items.push_back(itm); - } else { - heavy_items.push_back(itm); - } - } - - while (!light_items.empty() && !heavy_items.empty()) { - item l = light_items.back(); - item h = heavy_items.back(); - table_item tbl_itm = {l.weight, l.id, h.id}; - m_local_alias_table.push_back(tbl_itm); - h.weight = (h.weight + l.weight) - 1; - light_items.pop_back(); - if (h.weight < 1) { - light_items.push_back(h); - heavy_items.pop_back(); - } - } - - // Either heavy items or light_items is empty, need to flush the non empty - // vector and add them to the alias table with a p value of 1 - while (!heavy_items.empty()) { - item h = heavy_items.back(); - table_item tbl_itm = {1.0, h.id, 0}; - m_local_alias_table.push_back(tbl_itm); - heavy_items.pop_back(); - } - while (!light_items.empty()) { - item l = light_items.back(); - table_item tbl_itm = {1.0, l.id, 0}; - m_local_alias_table.push_back(tbl_itm); - light_items.pop_back(); - } - m_comm.barrier(); + items_to_send.clear(); + if (curr_weight != 0) { + items_to_send.push_back({local_item.id, curr_weight}); } + dest_rank++; + } else { + items_to_send.push_back(local_item); + curr_weight += local_item.weight; + } + } - public: - - template - void async_sample(Visitor visitor, const VisitorArgs &...args) { // Sample should be provided a lambda/functor that takes as an argument the item type - - auto sample_wrapper = [](auto ptr_a_tbl, const VisitorArgs &...args) { - std::uniform_int_distribution table_item_dist(0, ptr_a_tbl->m_local_alias_table.size()-1); - table_item tbl_itm = ptr_a_tbl->m_local_alias_table[table_item_dist(ptr_a_tbl->m_rng)]; - Item s; - if (tbl_itm.p == 1) { - s = tbl_itm.a; - } else { - std::uniform_real_distribution zero_one_dist(0.0, 1.0); - float f = zero_one_dist(ptr_a_tbl->m_rng); - if (f < tbl_itm.p) { - s = tbl_itm.a; - } else { - s = tbl_itm.b; - } - } - Visitor *vis = nullptr; - ygm::meta::apply_optional(*vis, std::make_tuple(ptr_a_tbl), std::forward_as_tuple(s, args...)); - }; - - uint32_t dest_rank = m_rank_dist(m_rng); - m_comm.async(dest_rank, sample_wrapper, pthis, std::forward(args)...); - } - - private: - ygm::comm& m_comm; - ygm::ygm_ptr pthis; - std::uniform_int_distribution m_rank_dist; - RNG& m_rng; - std::vector m_local_items; - std::vector m_local_alias_table; + // Need to handle items left in items to send + if (items_to_send.size() > 0 && curr_weight > 0.00001 && dest_rank < m_comm.size()) { // Account for floating point errors + m_comm.async(dest_rank, [](std::vector items, ygm_items_ptr new_items_ptr) { + new_items_ptr->insert(new_items_ptr->end(), items.begin(), items.end()); + }, items_to_send, ptr_new_items); + } + + m_comm.barrier(); + std::swap(new_local_items, m_local_items); + } + + bool is_balanced() { + // Not a perfect check. Will check if every rank has the same weight for each table + // but does not confirm that the actual target weight is met on each rank + double local_weight = 0.0; + for (const auto& itm : m_local_items) { + local_weight += itm.weight; + } + m_comm.barrier(); + auto equal = [](double a, double b){ + return (std::abs(a - b) < 0.00001); + }; + bool balanced = ygm::is_same(local_weight, m_comm, equal); + return balanced; + } + + void build_alias_table() { + // Make average weight of items 1, so coin flip is simple across all tables + double local_weight = 0.0; + for (const auto& itm : m_local_items) { + local_weight += itm.weight; + } + double avg_weight = local_weight / m_local_items.size(); + + double new_total_weight = 0; + for (auto& itm : m_local_items) { + itm.weight /= avg_weight; + new_total_weight += itm.weight; + } + double new_avg_weight = new_total_weight / m_local_items.size(); // Should be 1 + if (std::abs((new_avg_weight - 1.0)) < 0.0001) { + YGM_ASSERT_RELEASE(std::abs((new_avg_weight - 1.0)) < 0.00001); + } + + // Implementation of Vose's algorithm, utilized Keith Schwarz numerically stable version + // https://www.keithschwarz.com/darts-dice-coins/ + std::vector heavy_items; + std::vector light_items; + for (auto& itm : m_local_items) { + if (itm.weight < 1.0) { + light_items.push_back(itm); + } else { + heavy_items.push_back(itm); + } + } + + while (!light_items.empty() && !heavy_items.empty()) { + item l = light_items.back(); + item h = heavy_items.back(); + table_item tbl_itm = {l.weight, l.id, h.id}; + m_local_alias_table.push_back(tbl_itm); + h.weight = (h.weight + l.weight) - 1; + light_items.pop_back(); + if (h.weight < 1) { + light_items.push_back(h); + heavy_items.pop_back(); + } + } + + // Either heavy items or light_items is empty, need to flush the non empty + // vector and add them to the alias table with a p value of 1 + while (!heavy_items.empty()) { + item h = heavy_items.back(); + table_item tbl_itm = {1.0, h.id, 0}; + m_local_alias_table.push_back(tbl_itm); + heavy_items.pop_back(); + } + while (!light_items.empty()) { + item l = light_items.back(); + table_item tbl_itm = {1.0, l.id, 0}; + m_local_alias_table.push_back(tbl_itm); + light_items.pop_back(); + } + m_comm.barrier(); + } + + public: + + template + void async_sample([[maybe_unused]] Visitor visitor, const VisitorArgs &...args) { // Sample should be provided a lambda/functor that takes as an argument the item type + + auto sample_wrapper = [](auto ptr_a_tbl, const VisitorArgs &...args) { + std::uniform_int_distribution table_item_dist(0, ptr_a_tbl->m_local_alias_table.size()-1); + table_item tbl_itm = ptr_a_tbl->m_local_alias_table[table_item_dist(ptr_a_tbl->m_rng)]; + Item s; + if (tbl_itm.p == 1) { + s = tbl_itm.a; + } else { + // std::uniform_real_distribution zero_one_dist(0.0, 1.0); + float f = ptr_a_tbl->m_zero_one_dist(ptr_a_tbl->m_rng); + if (f < tbl_itm.p) { + s = tbl_itm.a; + } else { + s = tbl_itm.b; + } + } + Visitor *vis = nullptr; + ygm::meta::apply_optional(*vis, std::make_tuple(ptr_a_tbl), std::forward_as_tuple(s, args...)); + }; + + uint32_t dest_rank = m_rank_dist(m_rng); + m_comm.async(dest_rank, sample_wrapper, pthis, std::forward(args)...); + } + + private: + ygm::comm& m_comm; + ygm::ygm_ptr pthis; + RNG& m_rng; + std::uniform_int_distribution m_rank_dist; + std::vector m_local_items; + std::vector m_local_alias_table; + std::uniform_real_distribution m_zero_one_dist; }; } // namespace ygm::random \ No newline at end of file diff --git a/include/ygm/random/random.hpp b/include/ygm/random/random.hpp index 952acd54..6936c5fd 100644 --- a/include/ygm/random/random.hpp +++ b/include/ygm/random/random.hpp @@ -9,10 +9,49 @@ namespace ygm::random { + +/// @brief Applies a simple offset to the specified seed according to rank index +/// @tparam ResultType The random number (seed) type; defaults to std::mt19937 +/// @param comm The ygm::comm to be used +/// @param seed The specified seed +/// @return simply returns seed + rank +template +ResultType simple_offset(ygm::comm &comm, ResultType seed) { + return seed + comm.rank(); +} + +/// @brief A wrapper around a per-rank random engine that manipulates each +/// rank's seed according to a specified strategy +/// @tparam RandomEngine The underlying random engine, e.g. std::mt19337 +/// @tparam Function A `(ygm::comm, result_type) -> result_type` function that +/// modifies seeds for each rank +template +class random_engine { + public: + using rng_type = RandomEngine; + using result_type = typename RandomEngine::result_type; + + random_engine(ygm::comm &comm, result_type seed = std::random_device{}()) + : m_rng(Function(comm, seed)), m_seed(Function(comm, seed)) {} + + result_type operator()() { return m_rng(); } + + constexpr const result_type &seed() const { return m_seed; } + + static constexpr result_type min() { return rng_type::min(); } + static constexpr result_type max() { return rng_type::max(); } + + private: + rng_type m_rng; + result_type m_seed; +}; + /// @brief A simple offset rng alias /// @tparam RandomEngine The underlying random engine, e.g. std::mt19937 template using default_random_engine = - ygm::detail::random_engine; + random_engine; -} // namespace ygm \ No newline at end of file +} // namespace ygm::random \ No newline at end of file diff --git a/test/test_alias_table.cpp b/test/test_alias_table.cpp index 9617aff7..40188940 100644 --- a/test/test_alias_table.cpp +++ b/test/test_alias_table.cpp @@ -18,10 +18,10 @@ int main(int argc, char** argv) { int seed = 42; YGM_RNG ygm_rng(world, seed); + // // Test constructor works with YGM value container { ygm::container::bag> bag_of_items(world); - // ygm::container::bag> bag_of_items(world); uint32_t n_items_per_rank = 1000; const int max_item_weight = 100; @@ -37,11 +37,13 @@ int main(int argc, char** argv) { ygm::random::alias_table alias_tbl(world, ygm_rng, bag_of_items); - static uint32_t samples; + static uint32_t samples; + static uint32_t sum_of_samples; uint32_t samples_per_rank = 1000; - for (int i = 0; i < samples_per_rank; i++) { - alias_tbl.async_sample([](auto ptr, uint32_t item){ + for (uint32_t i = 0; i < samples_per_rank; i++) { + alias_tbl.async_sample([](uint32_t item){ samples++; + sum_of_samples += item; }); } world.barrier(); @@ -49,6 +51,7 @@ int main(int argc, char** argv) { YGM_ASSERT_RELEASE(total_samples == (samples_per_rank * world.size())); } + // // Test constructor works with YGM associative container { ygm::container::map map_of_items(world); @@ -69,7 +72,7 @@ int main(int argc, char** argv) { static uint32_t samples; uint32_t samples_per_rank = 1000; - for (int i = 0; i < samples_per_rank; i++) { + for (uint32_t i = 0; i < samples_per_rank; i++) { alias_tbl.async_sample([](auto ptr, uint32_t item){ samples++; }); From f571c5465af84f1289de268478f739438ac2a66c Mon Sep 17 00:00:00 2001 From: Lance Fletcher Date: Wed, 10 Dec 2025 15:48:26 -0800 Subject: [PATCH 05/15] Modified the traits of the alias_table constructor that takes in a YGM value container --- include/ygm/random/alias_table.hpp | 95 ++++++++++++++++-------------- test/test_alias_table.cpp | 54 ++++++++--------- 2 files changed, 78 insertions(+), 71 deletions(-) diff --git a/include/ygm/random/alias_table.hpp b/include/ygm/random/alias_table.hpp index 2113fbf1..9d27369e 100644 --- a/include/ygm/random/alias_table.hpp +++ b/include/ygm/random/alias_table.hpp @@ -7,20 +7,18 @@ #pragma once #include -#include #include #include -#include -#include #include #include -#include - namespace ygm::random { template -concept convertable_to_double = std::convertible_to; +concept second_within_convertable_to_double = + std::convertible_to>,double>; + +// convertable_to_double>> template concept is_associative_container = requires { @@ -28,17 +26,17 @@ concept is_associative_container = requires { typename T::mapped_type; }; -template -concept is_value_container = requires { - typename T::value_type; -}; - // template -// concept valid_stl_associative_container = requires { -// typename T::key_type; +// concept is_value_container = requires { // typename T::value_type; // }; +template +concept valid_stl_associative_container = requires { + typename T::key_type; + typename T::value_type; +}; + // template // concept valid_stl_value_container = requires { @@ -73,12 +71,20 @@ class alias_table { Item b; }; - template - requires ygm::container::detail::HasForAll && - is_value_container && - ygm::container::detail::DoubleItemTuple && - convertable_to_double> - alias_table(ygm::comm &comm, RNG &rng, Container &c) + // template + // requires ygm::container::detail::HasForAll && + // is_value_container && + // ygm::container::detail::DoubleItemTuple && + // convertable_to_double> + // template + // requires ygm::container::detail::DoubleItemTuple && + // convertable_to_double> + + template + requires ygm::container::detail::HasForAll && + ygm::container::detail::SingleItemTuple && + second_within_convertable_to_double + alias_table(ygm::comm &comm, RNG &rng, YGMContainer &c) : m_comm(comm), pthis(this), m_rng(rng), m_rank_dist(0, comm.size()-1), m_zero_one_dist(0.0, 1.0) { @@ -93,24 +99,24 @@ class alias_table { m_local_items.clear(); } - template - requires ygm::container::detail::HasForAll && - is_associative_container && - convertable_to_double - alias_table(ygm::comm &comm, RNG &rng, Container &c) - : m_comm(comm), pthis(this), m_rng(rng), - m_rank_dist(0, comm.size()-1), m_zero_one_dist(0.0, 1.0) { - - c.for_all([&](const auto &key, const auto &value){ - m_local_items.push_back({key, value}); - }); - comm.barrier(); - balance_weight(); - comm.barrier(); - YGM_ASSERT_RELEASE(is_balanced()); - build_alias_table(); - m_local_items.clear(); - } + // // template + // template + // requires ygm::container::detail::SingleItemTuple && + // convertable_to_double> + // alias_table(ygm::comm &comm, RNG &rng, YGMContainer &c) + // : m_comm(comm), pthis(this), m_rng(rng), + // m_rank_dist(0, comm.size()-1), m_zero_one_dist(0.0, 1.0) { + + // c.for_all([&](const auto &key, const auto &value){ + // m_local_items.push_back({key, value}); + // }); + // comm.barrier(); + // balance_weight(); + // comm.barrier(); + // YGM_ASSERT_RELEASE(is_balanced()); + // build_alias_table(); + // m_local_items.clear(); + // } private: @@ -172,8 +178,8 @@ class alias_table { } } - // Need to handle items left in items to send - if (items_to_send.size() > 0 && curr_weight > 0.00001 && dest_rank < m_comm.size()) { // Account for floating point errors + // Need to handle items left in items to send. Must also account for floating point errors. + if (items_to_send.size() > 0 && curr_weight > 0.00001 && dest_rank < m_comm.size()) { m_comm.async(dest_rank, [](std::vector items, ygm_items_ptr new_items_ptr) { new_items_ptr->insert(new_items_ptr->end(), items.begin(), items.end()); }, items_to_send, ptr_new_items); @@ -199,7 +205,7 @@ class alias_table { } void build_alias_table() { - // Make average weight of items 1, so coin flip is simple across all tables + // Make average weight of items 1. Makes random number generated for each sample be in [0,1] double local_weight = 0.0; for (const auto& itm : m_local_items) { local_weight += itm.weight; @@ -256,6 +262,7 @@ class alias_table { light_items.pop_back(); } m_comm.barrier(); + m_num_items_uniform_dist = std::uniform_int_distribution(0,m_local_alias_table.size()-1); } public: @@ -264,13 +271,11 @@ class alias_table { void async_sample([[maybe_unused]] Visitor visitor, const VisitorArgs &...args) { // Sample should be provided a lambda/functor that takes as an argument the item type auto sample_wrapper = [](auto ptr_a_tbl, const VisitorArgs &...args) { - std::uniform_int_distribution table_item_dist(0, ptr_a_tbl->m_local_alias_table.size()-1); - table_item tbl_itm = ptr_a_tbl->m_local_alias_table[table_item_dist(ptr_a_tbl->m_rng)]; + table_item tbl_itm = ptr_a_tbl->m_local_alias_table[ptr_a_tbl->m_num_items_uniform_dist(ptr_a_tbl->m_rng)]; Item s; if (tbl_itm.p == 1) { s = tbl_itm.a; } else { - // std::uniform_real_distribution zero_one_dist(0.0, 1.0); float f = ptr_a_tbl->m_zero_one_dist(ptr_a_tbl->m_rng); if (f < tbl_itm.p) { s = tbl_itm.a; @@ -290,9 +295,11 @@ class alias_table { ygm::comm& m_comm; ygm::ygm_ptr pthis; RNG& m_rng; - std::uniform_int_distribution m_rank_dist; std::vector m_local_items; std::vector m_local_alias_table; + + std::uniform_int_distribution m_rank_dist; + std::uniform_int_distribution m_num_items_uniform_dist; std::uniform_real_distribution m_zero_one_dist; }; diff --git a/test/test_alias_table.cpp b/test/test_alias_table.cpp index 40188940..cfcb12d8 100644 --- a/test/test_alias_table.cpp +++ b/test/test_alias_table.cpp @@ -19,7 +19,7 @@ int main(int argc, char** argv) { YGM_RNG ygm_rng(world, seed); // - // Test constructor works with YGM value container + // Test constructor works with YGM bag container { ygm::container::bag> bag_of_items(world); @@ -51,36 +51,36 @@ int main(int argc, char** argv) { YGM_ASSERT_RELEASE(total_samples == (samples_per_rank * world.size())); } - // - // Test constructor works with YGM associative container - { - ygm::container::map map_of_items(world); + // // + // // Test constructor works with YGM associative container + // { + // ygm::container::map map_of_items(world); - uint32_t n_items_per_rank = 1000; - const int max_item_weight = 100; + // uint32_t n_items_per_rank = 1000; + // const int max_item_weight = 100; - std::uniform_real_distribution dist(0, max_item_weight); - const uint32_t rank = world.rank(); - for (uint32_t i = 0; i < n_items_per_rank; i++) { - uint32_t id = i + rank * n_items_per_rank; - double w = dist(ygm_rng); - map_of_items.async_insert(id, w); - } - world.barrier(); + // std::uniform_real_distribution dist(0, max_item_weight); + // const uint32_t rank = world.rank(); + // for (uint32_t i = 0; i < n_items_per_rank; i++) { + // uint32_t id = i + rank * n_items_per_rank; + // double w = dist(ygm_rng); + // map_of_items.async_insert(id, w); + // } + // world.barrier(); - ygm::random::alias_table alias_tbl(world, ygm_rng, map_of_items); + // ygm::random::alias_table alias_tbl(world, ygm_rng, map_of_items); - static uint32_t samples; - uint32_t samples_per_rank = 1000; - for (uint32_t i = 0; i < samples_per_rank; i++) { - alias_tbl.async_sample([](auto ptr, uint32_t item){ - samples++; - }); - } - world.barrier(); - uint32_t total_samples = ygm::sum(samples, world); - YGM_ASSERT_RELEASE(total_samples == (samples_per_rank * world.size())); - } + // static uint32_t samples; + // uint32_t samples_per_rank = 1000; + // for (uint32_t i = 0; i < samples_per_rank; i++) { + // alias_tbl.async_sample([](auto ptr, uint32_t item){ + // samples++; + // }); + // } + // world.barrier(); + // uint32_t total_samples = ygm::sum(samples, world); + // YGM_ASSERT_RELEASE(total_samples == (samples_per_rank * world.size())); + // } return 0; } From 36b8d7eacd69bedd06af47e2a386366325856ec9 Mon Sep 17 00:00:00 2001 From: Lance Fletcher Date: Thu, 11 Dec 2025 08:46:11 -0800 Subject: [PATCH 06/15] Updated concepts to make use of existing YGM SingleItemTuple DoubleItemTuple concepts that take as input for_all_args type. Also added constructor that builds the alias table from an STL container --- include/ygm/random/alias_table.hpp | 119 +++++++++++++---------------- test/test_alias_table.cpp | 52 ++++++------- 2 files changed, 78 insertions(+), 93 deletions(-) diff --git a/include/ygm/random/alias_table.hpp b/include/ygm/random/alias_table.hpp index 9d27369e..90d9a5a7 100644 --- a/include/ygm/random/alias_table.hpp +++ b/include/ygm/random/alias_table.hpp @@ -15,42 +15,25 @@ namespace ygm::random { template -concept second_within_convertable_to_double = +concept second_within_convertible_to_double = std::convertible_to>,double>; -// convertable_to_double>> - -template -concept is_associative_container = requires { - typename T::key_type; - typename T::mapped_type; -}; - -// template -// concept is_value_container = requires { -// typename T::value_type; -// }; - -template -concept valid_stl_associative_container = requires { - typename T::key_type; - typename T::value_type; -}; - - -// template -// concept valid_stl_value_container = requires { -// typename T::value_type; -// requires ygm::container::detail::DoubleItemTuple; -// }; - -// template -// concept valid_stl_container = -// is_associative_container || valid_stl_value_container; +template +concept pair_like_and_convertible_to_weighted_item = + (requires { typename T::first_type; typename T::second_type; } && + std::is_convertible_v>) || + (std::tuple_size_v == 2 && + std::is_convertible_v, Item> && + std::is_convertible_v, double>); +// template +// concept second_within_convertible_to_double = +// std::convertible_to>,double>; template class alias_table { + + public: using self_type = alias_table; using weight_t = double; @@ -71,55 +54,57 @@ class alias_table { Item b; }; - // template - // requires ygm::container::detail::HasForAll && - // is_value_container && - // ygm::container::detail::DoubleItemTuple && - // convertable_to_double> - // template - // requires ygm::container::detail::DoubleItemTuple && - // convertable_to_double> - template requires ygm::container::detail::HasForAll && - ygm::container::detail::SingleItemTuple && - second_within_convertable_to_double + ygm::container::detail::SingleItemTuple && + pair_like_and_convertible_to_weighted_item> alias_table(ygm::comm &comm, RNG &rng, YGMContainer &c) : m_comm(comm), pthis(this), m_rng(rng), m_rank_dist(0, comm.size()-1), m_zero_one_dist(0.0, 1.0) { - c.for_all([&](const auto& id_weight_pair){ m_local_items.push_back({std::get<0>(id_weight_pair), std::get<1>(id_weight_pair)}); }); - comm.barrier(); - balance_weight(); - comm.barrier(); - YGM_ASSERT_RELEASE(is_balanced()); build_alias_table(); - m_local_items.clear(); } - // // template - // template - // requires ygm::container::detail::SingleItemTuple && - // convertable_to_double> - // alias_table(ygm::comm &comm, RNG &rng, YGMContainer &c) - // : m_comm(comm), pthis(this), m_rng(rng), - // m_rank_dist(0, comm.size()-1), m_zero_one_dist(0.0, 1.0) { - - // c.for_all([&](const auto &key, const auto &value){ - // m_local_items.push_back({key, value}); - // }); - // comm.barrier(); - // balance_weight(); - // comm.barrier(); - // YGM_ASSERT_RELEASE(is_balanced()); - // build_alias_table(); - // m_local_items.clear(); - // } + template + requires ygm::container::detail::HasForAll && + ygm::container::detail::DoubleItemTuple && + pair_like_and_convertible_to_weighted_item + alias_table(ygm::comm &comm, RNG &rng, YGMContainer &c) + : m_comm(comm), pthis(this), m_rng(rng), + m_rank_dist(0, comm.size()-1), m_zero_one_dist(0.0, 1.0) { + c.for_all([&](const auto &id, const auto &weight){ + m_local_items.push_back({id, weight}); + }); + build_alias_table(); + } + + template + requires ygm::container::detail::STLContainer && + pair_like_and_convertible_to_weighted_item< + Item, typename STLContainer::value_type> + alias_table(ygm::comm &comm, RNG &rng, STLContainer &c) + : m_comm(comm), pthis(this), m_rng(rng), + m_rank_dist(0, comm.size()-1), m_zero_one_dist(0.0, 1.0) { + for (const auto& [id, weight] : c) { + m_local_items.push_back({id, weight}); + } + build_alias_table(); + } private: + void build_alias_table() { + m_comm.barrier(); + balance_weight(); + m_comm.barrier(); + YGM_ASSERT_RELEASE(is_balanced()); + build_local_alias_table(); + m_local_items.clear(); + } + void balance_weight() { MPI_Comm comm = m_comm.get_mpi_comm(); double local_weight = 0.0; @@ -135,7 +120,7 @@ class alias_table { double target_weight = global_weight / m_comm.size(); // target weight per rank - uint32_t dest_rank = prfx_sum_weight / target_weight; + int dest_rank = prfx_sum_weight / target_weight; double curr_weight = std::fmod(prfx_sum_weight, target_weight); // Spillage weight std::vector new_local_items; @@ -204,7 +189,7 @@ class alias_table { return balanced; } - void build_alias_table() { + void build_local_alias_table() { // Make average weight of items 1. Makes random number generated for each sample be in [0,1] double local_weight = 0.0; for (const auto& itm : m_local_items) { diff --git a/test/test_alias_table.cpp b/test/test_alias_table.cpp index cfcb12d8..34c957db 100644 --- a/test/test_alias_table.cpp +++ b/test/test_alias_table.cpp @@ -51,36 +51,36 @@ int main(int argc, char** argv) { YGM_ASSERT_RELEASE(total_samples == (samples_per_rank * world.size())); } - // // - // // Test constructor works with YGM associative container - // { - // ygm::container::map map_of_items(world); + // + // Test constructor works with YGM associative container + { + ygm::container::map map_of_items(world); - // uint32_t n_items_per_rank = 1000; - // const int max_item_weight = 100; + uint32_t n_items_per_rank = 1000; + const int max_item_weight = 100; - // std::uniform_real_distribution dist(0, max_item_weight); - // const uint32_t rank = world.rank(); - // for (uint32_t i = 0; i < n_items_per_rank; i++) { - // uint32_t id = i + rank * n_items_per_rank; - // double w = dist(ygm_rng); - // map_of_items.async_insert(id, w); - // } - // world.barrier(); + std::uniform_real_distribution dist(0, max_item_weight); + const uint32_t rank = world.rank(); + for (uint32_t i = 0; i < n_items_per_rank; i++) { + uint32_t id = i + rank * n_items_per_rank; + double w = dist(ygm_rng); + map_of_items.async_insert(id, w); + } + world.barrier(); - // ygm::random::alias_table alias_tbl(world, ygm_rng, map_of_items); + ygm::random::alias_table alias_tbl(world, ygm_rng, map_of_items); - // static uint32_t samples; - // uint32_t samples_per_rank = 1000; - // for (uint32_t i = 0; i < samples_per_rank; i++) { - // alias_tbl.async_sample([](auto ptr, uint32_t item){ - // samples++; - // }); - // } - // world.barrier(); - // uint32_t total_samples = ygm::sum(samples, world); - // YGM_ASSERT_RELEASE(total_samples == (samples_per_rank * world.size())); - // } + static uint32_t samples; + uint32_t samples_per_rank = 1000; + for (uint32_t i = 0; i < samples_per_rank; i++) { + alias_tbl.async_sample([]([[maybe_unused]] auto ptr, [[maybe_unused]] uint32_t item){ + samples++; + }); + } + world.barrier(); + uint32_t total_samples = ygm::sum(samples, world); + YGM_ASSERT_RELEASE(total_samples == (samples_per_rank * world.size())); + } return 0; } From 08cff2cf7a4b87909af18b5dbe79b23084848d7f Mon Sep 17 00:00:00 2001 From: Lance Fletcher Date: Thu, 11 Dec 2025 09:44:11 -0800 Subject: [PATCH 07/15] Modified async_sample function to now capture visitor lambda in wrapper lambda. --- include/ygm/random/alias_table.hpp | 24 +++++++++--------------- 1 file changed, 9 insertions(+), 15 deletions(-) diff --git a/include/ygm/random/alias_table.hpp b/include/ygm/random/alias_table.hpp index 90d9a5a7..c4e7e3d2 100644 --- a/include/ygm/random/alias_table.hpp +++ b/include/ygm/random/alias_table.hpp @@ -32,15 +32,12 @@ concept pair_like_and_convertible_to_weighted_item = template class alias_table { - - public: using self_type = alias_table; - using weight_t = double; struct item { Item id; - weight_t weight; + double weight; template void serialize(Archive& ar) { // Needed for cereal serialization @@ -138,7 +135,7 @@ class alias_table { item item_to_send = {local_item.id, weight_to_send}; items_to_send.push_back(item_to_send); - if ((curr_weight > 0.0001) && (dest_rank < m_comm.size())) { // Accounts for rounding errors + if ((curr_weight > m_tolerance) && (dest_rank < m_comm.size())) { // Accounts for rounding errors // Moves weights to dest_rank's new_local_items m_comm.async(dest_rank, [](std::vector items, ygm_items_ptr new_items_ptr) { new_items_ptr->insert(new_items_ptr->end(), items.begin(), items.end()); @@ -164,7 +161,7 @@ class alias_table { } // Need to handle items left in items to send. Must also account for floating point errors. - if (items_to_send.size() > 0 && curr_weight > 0.00001 && dest_rank < m_comm.size()) { + if (items_to_send.size() > 0 && curr_weight > m_tolerance && dest_rank < m_comm.size()) { m_comm.async(dest_rank, [](std::vector items, ygm_items_ptr new_items_ptr) { new_items_ptr->insert(new_items_ptr->end(), items.begin(), items.end()); }, items_to_send, ptr_new_items); @@ -183,7 +180,7 @@ class alias_table { } m_comm.barrier(); auto equal = [](double a, double b){ - return (std::abs(a - b) < 0.00001); + return (std::abs(a - b) < m_tolerance); }; bool balanced = ygm::is_same(local_weight, m_comm, equal); return balanced; @@ -203,9 +200,7 @@ class alias_table { new_total_weight += itm.weight; } double new_avg_weight = new_total_weight / m_local_items.size(); // Should be 1 - if (std::abs((new_avg_weight - 1.0)) < 0.0001) { - YGM_ASSERT_RELEASE(std::abs((new_avg_weight - 1.0)) < 0.00001); - } + YGM_ASSERT_RELEASE(std::abs((new_avg_weight - 1.0)) < m_tolerance); // Implementation of Vose's algorithm, utilized Keith Schwarz numerically stable version // https://www.keithschwarz.com/darts-dice-coins/ @@ -253,9 +248,9 @@ class alias_table { public: template - void async_sample([[maybe_unused]] Visitor visitor, const VisitorArgs &...args) { // Sample should be provided a lambda/functor that takes as an argument the item type + void async_sample(Visitor&& visitor, const VisitorArgs &...args) { - auto sample_wrapper = [](auto ptr_a_tbl, const VisitorArgs &...args) { + auto sample_wrapper = [visitor](auto ptr_a_tbl, const VisitorArgs &...args) { table_item tbl_itm = ptr_a_tbl->m_local_alias_table[ptr_a_tbl->m_num_items_uniform_dist(ptr_a_tbl->m_rng)]; Item s; if (tbl_itm.p == 1) { @@ -268,8 +263,7 @@ class alias_table { s = tbl_itm.b; } } - Visitor *vis = nullptr; - ygm::meta::apply_optional(*vis, std::make_tuple(ptr_a_tbl), std::forward_as_tuple(s, args...)); + ygm::meta::apply_optional(visitor, std::make_tuple(ptr_a_tbl), std::forward_as_tuple(s, args...)); }; uint32_t dest_rank = m_rank_dist(m_rng); @@ -282,10 +276,10 @@ class alias_table { RNG& m_rng; std::vector m_local_items; std::vector m_local_alias_table; - std::uniform_int_distribution m_rank_dist; std::uniform_int_distribution m_num_items_uniform_dist; std::uniform_real_distribution m_zero_one_dist; + double m_tolerance = 1e-6; }; } // namespace ygm::random \ No newline at end of file From fc7fde5490db05e1c72e8a631169af60fd7ed936 Mon Sep 17 00:00:00 2001 From: Lance Fletcher Date: Thu, 11 Dec 2025 15:19:22 -0800 Subject: [PATCH 08/15] Fixed local alias table construction bug. Simplified async_sample() logic. Wrote test for various constructors and tested sampling accuracy via sampling words from corpus and comparing to ground truth frequency --- include/ygm/random/alias_table.hpp | 77 ++++++------- test/test_alias_table.cpp | 168 ++++++++++++++++++++++------- 2 files changed, 167 insertions(+), 78 deletions(-) diff --git a/include/ygm/random/alias_table.hpp b/include/ygm/random/alias_table.hpp index c4e7e3d2..d23ded16 100644 --- a/include/ygm/random/alias_table.hpp +++ b/include/ygm/random/alias_table.hpp @@ -14,10 +14,6 @@ namespace ygm::random { -template -concept second_within_convertible_to_double = - std::convertible_to>,double>; - template concept pair_like_and_convertible_to_weighted_item = (requires { typename T::first_type; typename T::second_type; } && @@ -25,9 +21,6 @@ concept pair_like_and_convertible_to_weighted_item = (std::tuple_size_v == 2 && std::is_convertible_v, Item> && std::is_convertible_v, double>); -// template -// concept second_within_convertible_to_double = -// std::convertible_to>,double>; template class alias_table { @@ -40,13 +33,14 @@ class alias_table { double weight; template - void serialize(Archive& ar) { // Needed for cereal serialization + void serialize(Archive& ar) { ar(id, weight); } }; struct table_item { - double p; // prob item a is selected = p, prob item b is selected = (1 - p) + // prob item a is selected = p, prob item b is selected = (1 - p) + double p; Item a; Item b; }; @@ -58,9 +52,10 @@ class alias_table { std::tuple_element_t<0,typename YGMContainer::for_all_args>> alias_table(ygm::comm &comm, RNG &rng, YGMContainer &c) : m_comm(comm), pthis(this), m_rng(rng), - m_rank_dist(0, comm.size()-1), m_zero_one_dist(0.0, 1.0) { + m_rank_dist(0, comm.size()-1), + m_zero_one_dist(0.0, 1.0), m_tolerance(1e-4) { c.for_all([&](const auto& id_weight_pair){ - m_local_items.push_back({std::get<0>(id_weight_pair), std::get<1>(id_weight_pair)}); + m_local_items.emplace_back(std::get<0>(id_weight_pair), std::get<1>(id_weight_pair)); }); build_alias_table(); } @@ -71,9 +66,10 @@ class alias_table { pair_like_and_convertible_to_weighted_item alias_table(ygm::comm &comm, RNG &rng, YGMContainer &c) : m_comm(comm), pthis(this), m_rng(rng), - m_rank_dist(0, comm.size()-1), m_zero_one_dist(0.0, 1.0) { + m_rank_dist(0, comm.size()-1), + m_zero_one_dist(0.0, 1.0), m_tolerance(1e-4) { c.for_all([&](const auto &id, const auto &weight){ - m_local_items.push_back({id, weight}); + m_local_items.emplace_back(id,weight); }); build_alias_table(); } @@ -84,9 +80,10 @@ class alias_table { Item, typename STLContainer::value_type> alias_table(ygm::comm &comm, RNG &rng, STLContainer &c) : m_comm(comm), pthis(this), m_rng(rng), - m_rank_dist(0, comm.size()-1), m_zero_one_dist(0.0, 1.0) { + m_rank_dist(0, comm.size()-1), + m_zero_one_dist(0.0, 1.0), m_tolerance(1e-4) { for (const auto& [id, weight] : c) { - m_local_items.push_back({id, weight}); + m_local_items.emplace_back(id, weight); } build_alias_table(); } @@ -108,25 +105,25 @@ class alias_table { for (uint32_t i = 0; i < m_local_items.size(); i++) { local_weight += m_local_items[i].weight; } - double global_weight; MPI_Allreduce(&local_weight, &global_weight, 1, MPI_DOUBLE, MPI_SUM, comm); - double prfx_sum_weight; MPI_Exscan(&local_weight, &prfx_sum_weight, 1, MPI_DOUBLE, MPI_SUM, comm); - - double target_weight = global_weight / m_comm.size(); // target weight per rank + // target_weight = Amount of weight each rank should have after balancing + double target_weight = global_weight / m_comm.size(); int dest_rank = prfx_sum_weight / target_weight; - double curr_weight = std::fmod(prfx_sum_weight, target_weight); // Spillage weight + // Spillage weight i.e. weight being contributed by other processors to dest's rank local distribution + double curr_weight = std::fmod(prfx_sum_weight, target_weight); std::vector new_local_items; using ygm_items_ptr = ygm::ygm_ptr>; ygm_items_ptr ptr_new_items = m_comm.make_ygm_ptr(new_local_items); - m_comm.barrier(); // ensure ygm_ptr is created on every rank + m_comm.barrier(); std::vector items_to_send; - for (uint64_t i = 0; i < m_local_items.size(); i++) { // WARNING: size of m_local_items can grow during loop + // WARNING: size of m_local_items can grow during loop. Do not use iterators or pointers in the loop. + for (uint64_t i = 0; i < m_local_items.size(); i++) { item local_item = m_local_items[i]; if (curr_weight + local_item.weight >= target_weight) { double remaining_weight = curr_weight + local_item.weight - target_weight; @@ -179,8 +176,8 @@ class alias_table { local_weight += itm.weight; } m_comm.barrier(); - auto equal = [](double a, double b){ - return (std::abs(a - b) < m_tolerance); + auto equal = [this](double a, double b){ + return (std::abs(a - b) < this->m_tolerance); }; bool balanced = ygm::is_same(local_weight, m_comm, equal); return balanced; @@ -194,7 +191,7 @@ class alias_table { } double avg_weight = local_weight / m_local_items.size(); - double new_total_weight = 0; + double new_total_weight = 0.0; for (auto& itm : m_local_items) { itm.weight /= avg_weight; new_total_weight += itm.weight; @@ -215,13 +212,13 @@ class alias_table { } while (!light_items.empty() && !heavy_items.empty()) { - item l = light_items.back(); - item h = heavy_items.back(); + item& l = light_items.back(); + item& h = heavy_items.back(); table_item tbl_itm = {l.weight, l.id, h.id}; m_local_alias_table.push_back(tbl_itm); h.weight = (h.weight + l.weight) - 1; light_items.pop_back(); - if (h.weight < 1) { + if (h.weight < 1.0) { light_items.push_back(h); heavy_items.pop_back(); } @@ -230,14 +227,14 @@ class alias_table { // Either heavy items or light_items is empty, need to flush the non empty // vector and add them to the alias table with a p value of 1 while (!heavy_items.empty()) { - item h = heavy_items.back(); - table_item tbl_itm = {1.0, h.id, 0}; + item& h = heavy_items.back(); + table_item tbl_itm = {1.0, h.id, Item()}; m_local_alias_table.push_back(tbl_itm); heavy_items.pop_back(); } while (!light_items.empty()) { - item l = light_items.back(); - table_item tbl_itm = {1.0, l.id, 0}; + item& l = light_items.back(); + table_item tbl_itm = {1.0, l.id, Item()}; m_local_alias_table.push_back(tbl_itm); light_items.pop_back(); } @@ -252,14 +249,10 @@ class alias_table { auto sample_wrapper = [visitor](auto ptr_a_tbl, const VisitorArgs &...args) { table_item tbl_itm = ptr_a_tbl->m_local_alias_table[ptr_a_tbl->m_num_items_uniform_dist(ptr_a_tbl->m_rng)]; - Item s; - if (tbl_itm.p == 1) { - s = tbl_itm.a; - } else { - float f = ptr_a_tbl->m_zero_one_dist(ptr_a_tbl->m_rng); - if (f < tbl_itm.p) { - s = tbl_itm.a; - } else { + Item s = tbl_itm.a; + if (tbl_itm.p < 1) { + double f = ptr_a_tbl->m_zero_one_dist(ptr_a_tbl->m_rng); + if (f > tbl_itm.p) { s = tbl_itm.b; } } @@ -278,8 +271,8 @@ class alias_table { std::vector m_local_alias_table; std::uniform_int_distribution m_rank_dist; std::uniform_int_distribution m_num_items_uniform_dist; - std::uniform_real_distribution m_zero_one_dist; - double m_tolerance = 1e-6; + std::uniform_real_distribution m_zero_one_dist; + double m_tolerance; }; } // namespace ygm::random \ No newline at end of file diff --git a/test/test_alias_table.cpp b/test/test_alias_table.cpp index 34c957db..4ea23eb1 100644 --- a/test/test_alias_table.cpp +++ b/test/test_alias_table.cpp @@ -8,8 +8,13 @@ #include #include #include +#include +#include #include #include +#include +#include +#include int main(int argc, char** argv) { @@ -19,62 +24,84 @@ int main(int argc, char** argv) { YGM_RNG ygm_rng(world, seed); // - // Test constructor works with YGM bag container + // Testing various constructors { - ygm::container::bag> bag_of_items(world); - - uint32_t n_items_per_rank = 1000; + uint32_t n_items_per_rank = 10000; const int max_item_weight = 100; - std::uniform_real_distribution dist(0, max_item_weight); - const uint32_t rank = world.rank(); - for (uint32_t i = 0; i < n_items_per_rank; i++) { - uint32_t id = i + rank * n_items_per_rank; + { // Constructing from ygm::container::bag + ygm::container::bag> bag_of_items(world); + for (uint32_t i = 0; i < n_items_per_rank; i++) { + uint32_t id = world.rank() + i * world.size(); double w = dist(ygm_rng); bag_of_items.async_insert({id, w}); + } + world.barrier(); + ygm::random::alias_table alias_tbl(world, ygm_rng, bag_of_items); + } + { // Constructing from ygm::container::map + ygm::container::map map_of_items(world); + for (uint32_t i = 0; i < n_items_per_rank; i++) { + uint32_t id = world.rank() + i * world.size(); + double w = dist(ygm_rng); + map_of_items.async_insert({id, w}); + } + world.barrier(); + ygm::random::alias_table alias_tbl(world, ygm_rng, map_of_items); + } + { // Constructing from ygm::container::array + ygm::container::array array_of_weights(world, n_items_per_rank*world.size()); + for (uint32_t i = 0; i < n_items_per_rank; i++) { + uint32_t id = world.rank() + i * world.size(); + double w = dist(ygm_rng); + array_of_weights.async_set(id,w); + } + world.barrier(); + ygm::random::alias_table alias_tbl(world, ygm_rng, array_of_weights); + } + { // Constructing from std::vector + std::vector> vec_of_items; + for (uint32_t i = 0; i < n_items_per_rank; i++) { + uint32_t id = world.rank() + i * world.size(); + double w = dist(ygm_rng); + vec_of_items.push_back({id,w}); + } + world.barrier(); + ygm::random::alias_table alias_tbl(world, ygm_rng, vec_of_items); + } + { // Constructing from std::map + std::map items_map; + for (uint32_t i = 0; i < n_items_per_rank; i++) { + uint32_t id = world.rank() + i * world.size(); + double w = dist(ygm_rng); + items_map[id] = w; + } + world.barrier(); + ygm::random::alias_table alias_tbl(world, ygm_rng, items_map); } - world.barrier(); - - ygm::random::alias_table alias_tbl(world, ygm_rng, bag_of_items); - - static uint32_t samples; - static uint32_t sum_of_samples; - uint32_t samples_per_rank = 1000; - for (uint32_t i = 0; i < samples_per_rank; i++) { - alias_tbl.async_sample([](uint32_t item){ - samples++; - sum_of_samples += item; - }); - } - world.barrier(); - uint32_t total_samples = ygm::sum(samples, world); - YGM_ASSERT_RELEASE(total_samples == (samples_per_rank * world.size())); } // - // Test constructor works with YGM associative container + // Test sampling numbers { ygm::container::map map_of_items(world); - uint32_t n_items_per_rank = 1000; + uint32_t n_items_per_rank = 10000; const int max_item_weight = 100; - std::uniform_real_distribution dist(0, max_item_weight); - const uint32_t rank = world.rank(); for (uint32_t i = 0; i < n_items_per_rank; i++) { - uint32_t id = i + rank * n_items_per_rank; - double w = dist(ygm_rng); - map_of_items.async_insert(id, w); + uint32_t id = world.rank() + i * world.size(); + double w = dist(ygm_rng); + map_of_items.async_insert(id,w); } world.barrier(); - ygm::random::alias_table alias_tbl(world, ygm_rng, map_of_items); - static uint32_t samples; - uint32_t samples_per_rank = 1000; + static uint32_t samples = 0; + uint32_t samples_per_rank = 100000; for (uint32_t i = 0; i < samples_per_rank; i++) { - alias_tbl.async_sample([]([[maybe_unused]] auto ptr, [[maybe_unused]] uint32_t item){ - samples++; + alias_tbl.async_sample([]([[maybe_unused]] auto ptr, [[maybe_unused]] uint32_t item){ + samples++; }); } world.barrier(); @@ -82,5 +109,74 @@ int main(int argc, char** argv) { YGM_ASSERT_RELEASE(total_samples == (samples_per_rank * world.size())); } + // + // Test sampling words with probability proportional to their frequency in a corpus + { + std::vector words; + std::ifstream file("data/loremipsum/loremipsum_0.txt"); + ygm::container::counting_set word_counts(world); + + static std::string ipsum = "ipsum"; + uint32_t ipsum_count = 0; + static std::string sit = "sit"; + uint32_t sit_count = 0; + uint32_t total_words = 0; + if (world.rank0()) { + std::string word; + while (file >> word) { + word_counts.async_insert(word); + ++total_words; + if (word == ipsum) { + ++ipsum_count; + } else if (word == sit) { + ++sit_count; + } + } + } + ygm::random::alias_table alias_tbl(world, ygm_rng, word_counts); + world.barrier(); + file.close(); + + static uint32_t samples = 0; + static uint32_t sampled_ipsums = 0; + static uint32_t sampled_sits = 0; + uint32_t samples_per_rank = 100000; + for (uint32_t i = 0; i < samples_per_rank; i++) { + alias_tbl.async_sample([](std::string word_sample){ + samples++; + if (word_sample == ipsum) { + ++sampled_ipsums; + } else if (word_sample == sit) { + ++sampled_sits; + } + }); + } + world.barrier(); + uint32_t total_samples = ygm::sum(samples, world); + uint32_t total_ipsums = ygm::sum(sampled_ipsums, world); + uint32_t total_sits = ygm::sum(sampled_sits, world); + + YGM_ASSERT_RELEASE(total_samples == (samples_per_rank * world.size())); + + if (world.rank() == 0) { + double ipsum_freq = double(ipsum_count) / total_words; + double sit_freq = double(sit_count) / total_words; + double ipsum_sample_freq = double(total_ipsums) / total_samples; + double sit_sample_freq = double(total_sits) / total_samples; + + world.cout0("\"ipsum\" actual frequency: ", ipsum_freq); + world.cout0("\"ipsum\" sample frequency: ", ipsum_sample_freq); + double dif = std::abs(ipsum_sample_freq - ipsum_freq); + world.cout0("\"ipsum\" frequency difference: ", dif); + YGM_ASSERT_RELEASE(dif < 1e-4); + + world.cout0("\"sit\" actual frequency: ", sit_freq); + world.cout0("\"sit\" sample frequency: ", sit_sample_freq); + dif = std::abs(sit_sample_freq - sit_freq); + world.cout0("\"sit\" frequency difference: ", dif); + YGM_ASSERT_RELEASE(dif < 1e-4); + } + } + return 0; } From 9587ee80bbc7bc2602251b532337774758b7a3aa Mon Sep 17 00:00:00 2001 From: Lance Fletcher Date: Thu, 11 Dec 2025 15:30:33 -0800 Subject: [PATCH 09/15] Added more sampling in alias table test to get better sampling frequency accuracy --- test/test_alias_table.cpp | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/test/test_alias_table.cpp b/test/test_alias_table.cpp index 4ea23eb1..812fd5fb 100644 --- a/test/test_alias_table.cpp +++ b/test/test_alias_table.cpp @@ -140,7 +140,7 @@ int main(int argc, char** argv) { static uint32_t samples = 0; static uint32_t sampled_ipsums = 0; static uint32_t sampled_sits = 0; - uint32_t samples_per_rank = 100000; + uint32_t samples_per_rank = 1000000; for (uint32_t i = 0; i < samples_per_rank; i++) { alias_tbl.async_sample([](std::string word_sample){ samples++; @@ -168,13 +168,13 @@ int main(int argc, char** argv) { world.cout0("\"ipsum\" sample frequency: ", ipsum_sample_freq); double dif = std::abs(ipsum_sample_freq - ipsum_freq); world.cout0("\"ipsum\" frequency difference: ", dif); - YGM_ASSERT_RELEASE(dif < 1e-4); + YGM_ASSERT_RELEASE(dif < 1e-3); world.cout0("\"sit\" actual frequency: ", sit_freq); world.cout0("\"sit\" sample frequency: ", sit_sample_freq); dif = std::abs(sit_sample_freq - sit_freq); world.cout0("\"sit\" frequency difference: ", dif); - YGM_ASSERT_RELEASE(dif < 1e-4); + YGM_ASSERT_RELEASE(dif < 1e-3); } } From 00c6f127ba6926483e15920f5debbee958bee038 Mon Sep 17 00:00:00 2001 From: Lance Fletcher Date: Fri, 12 Dec 2025 09:46:10 -0800 Subject: [PATCH 10/15] Wrote additional tests to ensure alias table weight balancing code can handle a variety of distributions --- include/ygm/random/alias_table.hpp | 4 +- test/test_alias_table.cpp | 66 +++++++++++++++++++++++++++--- 2 files changed, 63 insertions(+), 7 deletions(-) diff --git a/include/ygm/random/alias_table.hpp b/include/ygm/random/alias_table.hpp index d23ded16..3e21dd91 100644 --- a/include/ygm/random/alias_table.hpp +++ b/include/ygm/random/alias_table.hpp @@ -132,7 +132,7 @@ class alias_table { item item_to_send = {local_item.id, weight_to_send}; items_to_send.push_back(item_to_send); - if ((curr_weight > m_tolerance) && (dest_rank < m_comm.size())) { // Accounts for rounding errors + if ((curr_weight > 1e-6) && (dest_rank < m_comm.size())) { // Accounts for rounding errors // Moves weights to dest_rank's new_local_items m_comm.async(dest_rank, [](std::vector items, ygm_items_ptr new_items_ptr) { new_items_ptr->insert(new_items_ptr->end(), items.begin(), items.end()); @@ -158,7 +158,7 @@ class alias_table { } // Need to handle items left in items to send. Must also account for floating point errors. - if (items_to_send.size() > 0 && curr_weight > m_tolerance && dest_rank < m_comm.size()) { + if (items_to_send.size() > 0 && curr_weight > 1e-6 && dest_rank < m_comm.size()) { m_comm.async(dest_rank, [](std::vector items, ygm_items_ptr new_items_ptr) { new_items_ptr->insert(new_items_ptr->end(), items.begin(), items.end()); }, items_to_send, ptr_new_items); diff --git a/test/test_alias_table.cpp b/test/test_alias_table.cpp index 812fd5fb..8927c496 100644 --- a/test/test_alias_table.cpp +++ b/test/test_alias_table.cpp @@ -26,8 +26,8 @@ int main(int argc, char** argv) { // // Testing various constructors { - uint32_t n_items_per_rank = 10000; - const int max_item_weight = 100; + uint32_t n_items_per_rank = 100; + const int max_item_weight = 10; std::uniform_real_distribution dist(0, max_item_weight); { // Constructing from ygm::container::bag ygm::container::bag> bag_of_items(world); @@ -81,13 +81,69 @@ int main(int argc, char** argv) { } } + // + // Testing the construction of many distributions. Balancing capabilities being tested. + { + uint32_t alias_tables_to_construct = 1000; + uint32_t n_items_per_rank = 10000; + { // Testing uniform weight distribution + std::uniform_int_distribution max_item_weight_dist(50, 100); + for (uint32_t i = 0; i < alias_tables_to_construct; i++) { + ygm::container::map map_of_items(world); + uint32_t max_item_weight = max_item_weight_dist(ygm_rng); + std::uniform_real_distribution weight_dist(0, max_item_weight); + for (uint32_t j = 0; j < n_items_per_rank; j++) { + uint32_t id = world.rank() + j * world.size(); + double w = weight_dist(ygm_rng); + map_of_items.async_insert(id,w); + } + world.barrier(); + ygm::random::alias_table alias_tbl(world, ygm_rng, map_of_items); + } + } + { // Testing normal weight distribution + std::uniform_int_distribution mean_dist(50, 100); + std::uniform_int_distribution std_dev_dist(5, 20); + for (uint32_t i = 0; i < alias_tables_to_construct; i++) { + ygm::container::map map_of_items(world); + uint32_t mean = mean_dist(ygm_rng); + uint32_t std_dev = std_dev_dist(ygm_rng); + std::normal_distribution weight_dist(mean, std_dev); + for (uint32_t j = 0; j < n_items_per_rank; j++) { + uint32_t id = world.rank() + j * world.size(); + double w = weight_dist(ygm_rng); + map_of_items.async_insert(id,w); + } + world.barrier(); + ygm::random::alias_table alias_tbl(world, ygm_rng, map_of_items); + } + } + { // Testing gamma weight distribution + std::uniform_real_distribution alpha_dist(0.1, 10); + std::uniform_real_distribution theta_dist(10, 100); + for (uint32_t i = 0; i < alias_tables_to_construct; i++) { + ygm::container::map map_of_items(world); + double alpha = alpha_dist(ygm_rng); + double theta = theta_dist(ygm_rng); + std::gamma_distribution weight_dist(alpha, theta); + for (uint32_t j = 0; j < n_items_per_rank; j++) { + uint32_t id = world.rank() + j * world.size(); + double w = weight_dist(ygm_rng); + map_of_items.async_insert(id,w); + } + world.barrier(); + ygm::random::alias_table alias_tbl(world, ygm_rng, map_of_items); + } + } + } + // // Test sampling numbers { ygm::container::map map_of_items(world); - uint32_t n_items_per_rank = 10000; - const int max_item_weight = 100; + uint32_t n_items_per_rank = 1000; + uint32_t max_item_weight = 100; std::uniform_real_distribution dist(0, max_item_weight); for (uint32_t i = 0; i < n_items_per_rank; i++) { uint32_t id = world.rank() + i * world.size(); @@ -108,7 +164,7 @@ int main(int argc, char** argv) { uint32_t total_samples = ygm::sum(samples, world); YGM_ASSERT_RELEASE(total_samples == (samples_per_rank * world.size())); } - + // // Test sampling words with probability proportional to their frequency in a corpus { From 5af1458bde4d66d201c903c7da387481c1d8d745 Mon Sep 17 00:00:00 2001 From: Lance Fletcher Date: Mon, 15 Dec 2025 10:34:01 -0800 Subject: [PATCH 11/15] Altered alias table construction to not normalize weight. This means the average local weight is no longer 1. This was done to decrease the errors associated with inexact floating point operations --- include/ygm/random/alias_table.hpp | 66 ++++++++++++++++++------------ test/test_alias_table.cpp | 14 +++++-- 2 files changed, 50 insertions(+), 30 deletions(-) diff --git a/include/ygm/random/alias_table.hpp b/include/ygm/random/alias_table.hpp index 3e21dd91..3e1abdee 100644 --- a/include/ygm/random/alias_table.hpp +++ b/include/ygm/random/alias_table.hpp @@ -39,7 +39,7 @@ class alias_table { }; struct table_item { - // prob item a is selected = p, prob item b is selected = (1 - p) + // p / m_avg_weight = prob item a is selected. 1 - (p/m_avg_weight) is prob b is selected. double p; Item a; Item b; @@ -53,7 +53,7 @@ class alias_table { alias_table(ygm::comm &comm, RNG &rng, YGMContainer &c) : m_comm(comm), pthis(this), m_rng(rng), m_rank_dist(0, comm.size()-1), - m_zero_one_dist(0.0, 1.0), m_tolerance(1e-4) { + m_bucket_weight_dist(0.0, 1.0), m_tolerance(1e-4) { c.for_all([&](const auto& id_weight_pair){ m_local_items.emplace_back(std::get<0>(id_weight_pair), std::get<1>(id_weight_pair)); }); @@ -67,7 +67,7 @@ class alias_table { alias_table(ygm::comm &comm, RNG &rng, YGMContainer &c) : m_comm(comm), pthis(this), m_rng(rng), m_rank_dist(0, comm.size()-1), - m_zero_one_dist(0.0, 1.0), m_tolerance(1e-4) { + m_bucket_weight_dist(0.0, 1.0), m_tolerance(1e-4) { c.for_all([&](const auto &id, const auto &weight){ m_local_items.emplace_back(id,weight); }); @@ -81,7 +81,7 @@ class alias_table { alias_table(ygm::comm &comm, RNG &rng, STLContainer &c) : m_comm(comm), pthis(this), m_rng(rng), m_rank_dist(0, comm.size()-1), - m_zero_one_dist(0.0, 1.0), m_tolerance(1e-4) { + m_bucket_weight_dist(0.0, 1.0), m_tolerance(1e-4) { for (const auto& [id, weight] : c) { m_local_items.emplace_back(id, weight); } @@ -94,7 +94,6 @@ class alias_table { m_comm.barrier(); balance_weight(); m_comm.barrier(); - YGM_ASSERT_RELEASE(is_balanced()); build_local_alias_table(); m_local_items.clear(); } @@ -166,45 +165,56 @@ class alias_table { m_comm.barrier(); std::swap(new_local_items, m_local_items); + + YGM_ASSERT_RELEASE(is_balanced(target_weight)); } - bool is_balanced() { - // Not a perfect check. Will check if every rank has the same weight for each table - // but does not confirm that the actual target weight is met on each rank + bool is_balanced(double target) { + + double local_weight = 0.0; for (const auto& itm : m_local_items) { local_weight += itm.weight; } + + double dif = std::abs(target - local_weight); + if (dif > this->m_tolerance) { + std::cerr << "rank " << m_comm.rank() << "***********local weight: " << local_weight << ", target: " << target << ", dif: " << dif << "*************\n"; + } + m_comm.barrier(); auto equal = [this](double a, double b){ - return (std::abs(a - b) < this->m_tolerance); + // return (std::abs(a - b) < this->m_tolerance); + double dif = std::abs(a - b); + if (dif > this->m_tolerance) { + std::cerr << "rank " << this->m_comm.rank() << "=================a: " << a << ", b: " << b << ", dif: " << dif << "==================\n"; + return false; + } else { + return true; + } + // return (std::abs(a - b) < this->m_tolerance); }; bool balanced = ygm::is_same(local_weight, m_comm, equal); + if (!balanced) { + std::cerr << "Rank: " << m_comm.rank() << ", local weight: " << local_weight << "\n"; + std::cerr << "Rank: " << m_comm.rank() << ", total local items: " << m_local_items.size() << "\n"; + } return balanced; } void build_local_alias_table() { - // Make average weight of items 1. Makes random number generated for each sample be in [0,1] double local_weight = 0.0; for (const auto& itm : m_local_items) { local_weight += itm.weight; } double avg_weight = local_weight / m_local_items.size(); - double new_total_weight = 0.0; - for (auto& itm : m_local_items) { - itm.weight /= avg_weight; - new_total_weight += itm.weight; - } - double new_avg_weight = new_total_weight / m_local_items.size(); // Should be 1 - YGM_ASSERT_RELEASE(std::abs((new_avg_weight - 1.0)) < m_tolerance); - // Implementation of Vose's algorithm, utilized Keith Schwarz numerically stable version // https://www.keithschwarz.com/darts-dice-coins/ std::vector heavy_items; std::vector light_items; for (auto& itm : m_local_items) { - if (itm.weight < 1.0) { + if (itm.weight < avg_weight) { light_items.push_back(itm); } else { heavy_items.push_back(itm); @@ -216,30 +226,32 @@ class alias_table { item& h = heavy_items.back(); table_item tbl_itm = {l.weight, l.id, h.id}; m_local_alias_table.push_back(tbl_itm); - h.weight = (h.weight + l.weight) - 1; + h.weight = (h.weight + l.weight) - avg_weight; light_items.pop_back(); - if (h.weight < 1.0) { + if (h.weight < avg_weight) { light_items.push_back(h); heavy_items.pop_back(); } } // Either heavy items or light_items is empty, need to flush the non empty - // vector and add them to the alias table with a p value of 1 + // vector and add them to the alias table with a p value of avg_weight while (!heavy_items.empty()) { item& h = heavy_items.back(); - table_item tbl_itm = {1.0, h.id, Item()}; + table_item tbl_itm = {avg_weight, h.id, Item()}; m_local_alias_table.push_back(tbl_itm); heavy_items.pop_back(); } while (!light_items.empty()) { item& l = light_items.back(); - table_item tbl_itm = {1.0, l.id, Item()}; + table_item tbl_itm = {avg_weight, l.id, Item()}; m_local_alias_table.push_back(tbl_itm); light_items.pop_back(); } m_comm.barrier(); m_num_items_uniform_dist = std::uniform_int_distribution(0,m_local_alias_table.size()-1); + m_bucket_weight_dist = std::uniform_real_distribution(0,avg_weight); + m_avg_weight = avg_weight; } public: @@ -250,8 +262,8 @@ class alias_table { auto sample_wrapper = [visitor](auto ptr_a_tbl, const VisitorArgs &...args) { table_item tbl_itm = ptr_a_tbl->m_local_alias_table[ptr_a_tbl->m_num_items_uniform_dist(ptr_a_tbl->m_rng)]; Item s = tbl_itm.a; - if (tbl_itm.p < 1) { - double f = ptr_a_tbl->m_zero_one_dist(ptr_a_tbl->m_rng); + if (tbl_itm.p < ptr_a_tbl->m_avg_weight) { + double f = ptr_a_tbl->m_bucket_weight_dist(ptr_a_tbl->m_rng); if (f > tbl_itm.p) { s = tbl_itm.b; } @@ -272,7 +284,9 @@ class alias_table { std::uniform_int_distribution m_rank_dist; std::uniform_int_distribution m_num_items_uniform_dist; std::uniform_real_distribution m_zero_one_dist; + std::uniform_real_distribution m_bucket_weight_dist; double m_tolerance; + double m_avg_weight; }; } // namespace ygm::random \ No newline at end of file diff --git a/test/test_alias_table.cpp b/test/test_alias_table.cpp index 8927c496..04e78199 100644 --- a/test/test_alias_table.cpp +++ b/test/test_alias_table.cpp @@ -26,8 +26,8 @@ int main(int argc, char** argv) { // // Testing various constructors { - uint32_t n_items_per_rank = 100; - const int max_item_weight = 10; + uint32_t n_items_per_rank = 1000; + const int max_item_weight = 100; std::uniform_real_distribution dist(0, max_item_weight); { // Constructing from ygm::container::bag ygm::container::bag> bag_of_items(world); @@ -84,11 +84,12 @@ int main(int argc, char** argv) { // // Testing the construction of many distributions. Balancing capabilities being tested. { - uint32_t alias_tables_to_construct = 1000; - uint32_t n_items_per_rank = 10000; + uint32_t alias_tables_to_construct = 20000; + uint32_t n_items_per_rank = 5000; { // Testing uniform weight distribution std::uniform_int_distribution max_item_weight_dist(50, 100); for (uint32_t i = 0; i < alias_tables_to_construct; i++) { + world.cout0("uniform distribution alias table ", i, " of ", alias_tables_to_construct); ygm::container::map map_of_items(world); uint32_t max_item_weight = max_item_weight_dist(ygm_rng); std::uniform_real_distribution weight_dist(0, max_item_weight); @@ -100,11 +101,13 @@ int main(int argc, char** argv) { world.barrier(); ygm::random::alias_table alias_tbl(world, ygm_rng, map_of_items); } + world.cout0("Finished uniform distribution alias table test"); } { // Testing normal weight distribution std::uniform_int_distribution mean_dist(50, 100); std::uniform_int_distribution std_dev_dist(5, 20); for (uint32_t i = 0; i < alias_tables_to_construct; i++) { + world.cout0("Normal distribution alias table ", i, " of ", alias_tables_to_construct); ygm::container::map map_of_items(world); uint32_t mean = mean_dist(ygm_rng); uint32_t std_dev = std_dev_dist(ygm_rng); @@ -117,11 +120,13 @@ int main(int argc, char** argv) { world.barrier(); ygm::random::alias_table alias_tbl(world, ygm_rng, map_of_items); } + world.cout0("Finished normal distribution alias table test"); } { // Testing gamma weight distribution std::uniform_real_distribution alpha_dist(0.1, 10); std::uniform_real_distribution theta_dist(10, 100); for (uint32_t i = 0; i < alias_tables_to_construct; i++) { + world.cout0("Gamma distribution alias table ", i, " of ", alias_tables_to_construct); ygm::container::map map_of_items(world); double alpha = alpha_dist(ygm_rng); double theta = theta_dist(ygm_rng); @@ -134,6 +139,7 @@ int main(int argc, char** argv) { world.barrier(); ygm::random::alias_table alias_tbl(world, ygm_rng, map_of_items); } + world.cout0("Finished gamma distribution alias table test"); } } From 5d1a065090c0e2bef6fd67442faf0bee50b300bd Mon Sep 17 00:00:00 2001 From: Lance Fletcher Date: Mon, 15 Dec 2025 12:55:32 -0800 Subject: [PATCH 12/15] Fixed uninitialized variable error in balancing function. Testing has showed this has fixed the bug, but consistently reproducing the bug has proved to be difficult, so bug might still be present --- include/ygm/random/alias_table.hpp | 34 +++++++----------------------- test/test_alias_table.cpp | 4 ++-- 2 files changed, 10 insertions(+), 28 deletions(-) diff --git a/include/ygm/random/alias_table.hpp b/include/ygm/random/alias_table.hpp index 3e1abdee..b3c8f475 100644 --- a/include/ygm/random/alias_table.hpp +++ b/include/ygm/random/alias_table.hpp @@ -104,9 +104,9 @@ class alias_table { for (uint32_t i = 0; i < m_local_items.size(); i++) { local_weight += m_local_items[i].weight; } - double global_weight; + double global_weight = 0; MPI_Allreduce(&local_weight, &global_weight, 1, MPI_DOUBLE, MPI_SUM, comm); - double prfx_sum_weight; + double prfx_sum_weight = 0; MPI_Exscan(&local_weight, &prfx_sum_weight, 1, MPI_DOUBLE, MPI_SUM, comm); // target_weight = Amount of weight each rank should have after balancing @@ -131,7 +131,7 @@ class alias_table { item item_to_send = {local_item.id, weight_to_send}; items_to_send.push_back(item_to_send); - if ((curr_weight > 1e-6) && (dest_rank < m_comm.size())) { // Accounts for rounding errors + if ((curr_weight > 1e-8) && (dest_rank < m_comm.size())) { // Accounts for rounding errors // Moves weights to dest_rank's new_local_items m_comm.async(dest_rank, [](std::vector items, ygm_items_ptr new_items_ptr) { new_items_ptr->insert(new_items_ptr->end(), items.begin(), items.end()); @@ -157,7 +157,7 @@ class alias_table { } // Need to handle items left in items to send. Must also account for floating point errors. - if (items_to_send.size() > 0 && curr_weight > 1e-6 && dest_rank < m_comm.size()) { + if (items_to_send.size() > 0 && curr_weight > 1e-8 && dest_rank < m_comm.size()) { m_comm.async(dest_rank, [](std::vector items, ygm_items_ptr new_items_ptr) { new_items_ptr->insert(new_items_ptr->end(), items.begin(), items.end()); }, items_to_send, ptr_new_items); @@ -166,39 +166,21 @@ class alias_table { m_comm.barrier(); std::swap(new_local_items, m_local_items); + YGM_ASSERT_RELEASE(m_local_items.size() > 0); YGM_ASSERT_RELEASE(is_balanced(target_weight)); } - bool is_balanced(double target) { - - + bool is_balanced(double target) { double local_weight = 0.0; for (const auto& itm : m_local_items) { local_weight += itm.weight; } - - double dif = std::abs(target - local_weight); - if (dif > this->m_tolerance) { - std::cerr << "rank " << m_comm.rank() << "***********local weight: " << local_weight << ", target: " << target << ", dif: " << dif << "*************\n"; - } - + YGM_ASSERT_RELEASE(std::abs(target - local_weight) < m_tolerance); m_comm.barrier(); auto equal = [this](double a, double b){ - // return (std::abs(a - b) < this->m_tolerance); - double dif = std::abs(a - b); - if (dif > this->m_tolerance) { - std::cerr << "rank " << this->m_comm.rank() << "=================a: " << a << ", b: " << b << ", dif: " << dif << "==================\n"; - return false; - } else { - return true; - } - // return (std::abs(a - b) < this->m_tolerance); + return (std::abs(a - b) < this->m_tolerance); }; bool balanced = ygm::is_same(local_weight, m_comm, equal); - if (!balanced) { - std::cerr << "Rank: " << m_comm.rank() << ", local weight: " << local_weight << "\n"; - std::cerr << "Rank: " << m_comm.rank() << ", total local items: " << m_local_items.size() << "\n"; - } return balanced; } diff --git a/test/test_alias_table.cpp b/test/test_alias_table.cpp index 04e78199..c2f1b850 100644 --- a/test/test_alias_table.cpp +++ b/test/test_alias_table.cpp @@ -84,8 +84,8 @@ int main(int argc, char** argv) { // // Testing the construction of many distributions. Balancing capabilities being tested. { - uint32_t alias_tables_to_construct = 20000; - uint32_t n_items_per_rank = 5000; + uint32_t alias_tables_to_construct = 10000; + uint32_t n_items_per_rank = 1000; { // Testing uniform weight distribution std::uniform_int_distribution max_item_weight_dist(50, 100); for (uint32_t i = 0; i < alias_tables_to_construct; i++) { From 2e4eb771351fcf8f844ef4353dca7b2babd65131 Mon Sep 17 00:00:00 2001 From: Lance Fletcher Date: Tue, 30 Dec 2025 14:04:17 -0800 Subject: [PATCH 13/15] Removed unnecessary weight checks, also altered constructor to either take a seed for the rng or create an rng with a rag with a random seed --- include/ygm/random/alias_table.hpp | 52 ++++++++++++++++++------------ test/test_alias_table.cpp | 39 +++++++++++----------- 2 files changed, 50 insertions(+), 41 deletions(-) diff --git a/include/ygm/random/alias_table.hpp b/include/ygm/random/alias_table.hpp index b3c8f475..c9de59e5 100644 --- a/include/ygm/random/alias_table.hpp +++ b/include/ygm/random/alias_table.hpp @@ -22,11 +22,11 @@ concept pair_like_and_convertible_to_weighted_item = std::is_convertible_v, Item> && std::is_convertible_v, double>); -template +template class alias_table { public: - using self_type = alias_table; + using self_type = alias_table; struct item { Item id; @@ -39,7 +39,7 @@ class alias_table { }; struct table_item { - // p / m_avg_weight = prob item a is selected. 1 - (p/m_avg_weight) is prob b is selected. + // p / m_avg_weight = prob item a is selected. 1 - p / m_avg_weight is prob b is selected. double p; Item a; Item b; @@ -50,10 +50,13 @@ class alias_table { ygm::container::detail::SingleItemTuple && pair_like_and_convertible_to_weighted_item> - alias_table(ygm::comm &comm, RNG &rng, YGMContainer &c) - : m_comm(comm), pthis(this), m_rng(rng), - m_rank_dist(0, comm.size()-1), - m_bucket_weight_dist(0.0, 1.0), m_tolerance(1e-4) { + alias_table(ygm::comm &comm, YGMContainer &c, + std::optional seed = std::nullopt) + : m_comm(comm), pthis(this), m_rank_dist(0, comm.size()-1), + m_bucket_weight_dist(0.0, 1.0), m_rng(comm) { + if (seed) { + m_rng = default_random_engine<>(comm, *seed); + } c.for_all([&](const auto& id_weight_pair){ m_local_items.emplace_back(std::get<0>(id_weight_pair), std::get<1>(id_weight_pair)); }); @@ -64,10 +67,13 @@ class alias_table { requires ygm::container::detail::HasForAll && ygm::container::detail::DoubleItemTuple && pair_like_and_convertible_to_weighted_item - alias_table(ygm::comm &comm, RNG &rng, YGMContainer &c) - : m_comm(comm), pthis(this), m_rng(rng), - m_rank_dist(0, comm.size()-1), - m_bucket_weight_dist(0.0, 1.0), m_tolerance(1e-4) { + alias_table(ygm::comm &comm, YGMContainer &c, + std::optional seed = std::nullopt) + : m_comm(comm), pthis(this), m_rank_dist(0, comm.size()-1), + m_bucket_weight_dist(0.0, 1.0), m_rng(comm) { + if (seed) { + m_rng = default_random_engine<>(comm, *seed); + } c.for_all([&](const auto &id, const auto &weight){ m_local_items.emplace_back(id,weight); }); @@ -78,10 +84,13 @@ class alias_table { requires ygm::container::detail::STLContainer && pair_like_and_convertible_to_weighted_item< Item, typename STLContainer::value_type> - alias_table(ygm::comm &comm, RNG &rng, STLContainer &c) - : m_comm(comm), pthis(this), m_rng(rng), - m_rank_dist(0, comm.size()-1), - m_bucket_weight_dist(0.0, 1.0), m_tolerance(1e-4) { + alias_table(ygm::comm &comm, STLContainer &c, + std::optional seed = std::nullopt) + : m_comm(comm), pthis(this), m_rank_dist(0, comm.size()-1), + m_bucket_weight_dist(0.0, 1.0), m_rng(comm) { + if (seed) { + m_rng = default_random_engine<>(comm, *seed); + } for (const auto& [id, weight] : c) { m_local_items.emplace_back(id, weight); } @@ -131,7 +140,7 @@ class alias_table { item item_to_send = {local_item.id, weight_to_send}; items_to_send.push_back(item_to_send); - if ((curr_weight > 1e-8) && (dest_rank < m_comm.size())) { // Accounts for rounding errors + if (dest_rank < m_comm.size()) { // Moves weights to dest_rank's new_local_items m_comm.async(dest_rank, [](std::vector items, ygm_items_ptr new_items_ptr) { new_items_ptr->insert(new_items_ptr->end(), items.begin(), items.end()); @@ -157,7 +166,7 @@ class alias_table { } // Need to handle items left in items to send. Must also account for floating point errors. - if (items_to_send.size() > 0 && curr_weight > 1e-8 && dest_rank < m_comm.size()) { + if (items_to_send.size() > 0 && dest_rank < m_comm.size()) { m_comm.async(dest_rank, [](std::vector items, ygm_items_ptr new_items_ptr) { new_items_ptr->insert(new_items_ptr->end(), items.begin(), items.end()); }, items_to_send, ptr_new_items); @@ -175,10 +184,12 @@ class alias_table { for (const auto& itm : m_local_items) { local_weight += itm.weight; } - YGM_ASSERT_RELEASE(std::abs(target - local_weight) < m_tolerance); + double dif = std::abs(target - local_weight); + YGM_ASSERT_RELEASE(dif < 1e-6); + m_comm.barrier(); auto equal = [this](double a, double b){ - return (std::abs(a - b) < this->m_tolerance); + return (std::abs(a - b) < 1e-6); }; bool balanced = ygm::is_same(local_weight, m_comm, equal); return balanced; @@ -260,15 +271,14 @@ class alias_table { private: ygm::comm& m_comm; ygm::ygm_ptr pthis; - RNG& m_rng; std::vector m_local_items; std::vector m_local_alias_table; std::uniform_int_distribution m_rank_dist; std::uniform_int_distribution m_num_items_uniform_dist; std::uniform_real_distribution m_zero_one_dist; std::uniform_real_distribution m_bucket_weight_dist; - double m_tolerance; double m_avg_weight; + ygm::random::default_random_engine<> m_rng; }; } // namespace ygm::random \ No newline at end of file diff --git a/test/test_alias_table.cpp b/test/test_alias_table.cpp index c2f1b850..1b0b8a51 100644 --- a/test/test_alias_table.cpp +++ b/test/test_alias_table.cpp @@ -19,9 +19,8 @@ int main(int argc, char** argv) { ygm::comm world(&argc, &argv); - using YGM_RNG = ygm::random::default_random_engine<>; int seed = 42; - YGM_RNG ygm_rng(world, seed); + ygm::random::default_random_engine<> ygm_rng(world, seed); // // Testing various constructors @@ -37,7 +36,7 @@ int main(int argc, char** argv) { bag_of_items.async_insert({id, w}); } world.barrier(); - ygm::random::alias_table alias_tbl(world, ygm_rng, bag_of_items); + ygm::random::alias_table alias_tbl(world, bag_of_items, ygm_rng()); } { // Constructing from ygm::container::map ygm::container::map map_of_items(world); @@ -47,7 +46,7 @@ int main(int argc, char** argv) { map_of_items.async_insert({id, w}); } world.barrier(); - ygm::random::alias_table alias_tbl(world, ygm_rng, map_of_items); + ygm::random::alias_table alias_tbl(world, map_of_items, ygm_rng()); } { // Constructing from ygm::container::array ygm::container::array array_of_weights(world, n_items_per_rank*world.size()); @@ -57,7 +56,7 @@ int main(int argc, char** argv) { array_of_weights.async_set(id,w); } world.barrier(); - ygm::random::alias_table alias_tbl(world, ygm_rng, array_of_weights); + ygm::random::alias_table alias_tbl(world, array_of_weights, ygm_rng()); } { // Constructing from std::vector std::vector> vec_of_items; @@ -67,7 +66,7 @@ int main(int argc, char** argv) { vec_of_items.push_back({id,w}); } world.barrier(); - ygm::random::alias_table alias_tbl(world, ygm_rng, vec_of_items); + ygm::random::alias_table alias_tbl(world, vec_of_items, ygm_rng()); } { // Constructing from std::map std::map items_map; @@ -77,14 +76,14 @@ int main(int argc, char** argv) { items_map[id] = w; } world.barrier(); - ygm::random::alias_table alias_tbl(world, ygm_rng, items_map); + ygm::random::alias_table alias_tbl(world, items_map, ygm_rng()); } } // // Testing the construction of many distributions. Balancing capabilities being tested. { - uint32_t alias_tables_to_construct = 10000; + uint32_t alias_tables_to_construct = 1000; uint32_t n_items_per_rank = 1000; { // Testing uniform weight distribution std::uniform_int_distribution max_item_weight_dist(50, 100); @@ -99,7 +98,7 @@ int main(int argc, char** argv) { map_of_items.async_insert(id,w); } world.barrier(); - ygm::random::alias_table alias_tbl(world, ygm_rng, map_of_items); + ygm::random::alias_table alias_tbl(world, map_of_items); } world.cout0("Finished uniform distribution alias table test"); } @@ -118,7 +117,7 @@ int main(int argc, char** argv) { map_of_items.async_insert(id,w); } world.barrier(); - ygm::random::alias_table alias_tbl(world, ygm_rng, map_of_items); + ygm::random::alias_table alias_tbl(world, map_of_items, ygm_rng()); } world.cout0("Finished normal distribution alias table test"); } @@ -137,7 +136,7 @@ int main(int argc, char** argv) { map_of_items.async_insert(id,w); } world.barrier(); - ygm::random::alias_table alias_tbl(world, ygm_rng, map_of_items); + ygm::random::alias_table alias_tbl(world, map_of_items, ygm_rng()); } world.cout0("Finished gamma distribution alias table test"); } @@ -157,7 +156,7 @@ int main(int argc, char** argv) { map_of_items.async_insert(id,w); } world.barrier(); - ygm::random::alias_table alias_tbl(world, ygm_rng, map_of_items); + ygm::random::alias_table alias_tbl(world, map_of_items, ygm_rng()); static uint32_t samples = 0; uint32_t samples_per_rank = 100000; @@ -195,14 +194,14 @@ int main(int argc, char** argv) { } } } - ygm::random::alias_table alias_tbl(world, ygm_rng, word_counts); + ygm::random::alias_table alias_tbl(world, word_counts); world.barrier(); file.close(); - static uint32_t samples = 0; - static uint32_t sampled_ipsums = 0; - static uint32_t sampled_sits = 0; - uint32_t samples_per_rank = 1000000; + static uint64_t samples = 0; + static uint64_t sampled_ipsums = 0; + static uint64_t sampled_sits = 0; + uint32_t samples_per_rank = 10000000; for (uint32_t i = 0; i < samples_per_rank; i++) { alias_tbl.async_sample([](std::string word_sample){ samples++; @@ -214,9 +213,9 @@ int main(int argc, char** argv) { }); } world.barrier(); - uint32_t total_samples = ygm::sum(samples, world); - uint32_t total_ipsums = ygm::sum(sampled_ipsums, world); - uint32_t total_sits = ygm::sum(sampled_sits, world); + uint64_t total_samples = ygm::sum(samples, world); + uint64_t total_ipsums = ygm::sum(sampled_ipsums, world); + uint64_t total_sits = ygm::sum(sampled_sits, world); YGM_ASSERT_RELEASE(total_samples == (samples_per_rank * world.size())); From 1ae9aa12495809b1d90d6895f1253b759a724485 Mon Sep 17 00:00:00 2001 From: Lance Fletcher Date: Tue, 30 Dec 2025 15:45:50 -0800 Subject: [PATCH 14/15] Modified test to adhere to new ygm::random organization --- test/test_multi_compilation_unit_lib.h | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/test/test_multi_compilation_unit_lib.h b/test/test_multi_compilation_unit_lib.h index 9af744d4..fef29903 100644 --- a/test/test_multi_compilation_unit_lib.h +++ b/test/test_multi_compilation_unit_lib.h @@ -1,7 +1,9 @@ // Including all of YGM (except parquet_parser) #include -#include + +#include +#include #include #include From 2557d4981f064c7724d1166b2a799ac7923cfbd9 Mon Sep 17 00:00:00 2001 From: Lance Fletcher Date: Mon, 5 Jan 2026 11:29:36 -0800 Subject: [PATCH 15/15] Changed implemenation to use ygm collectives rather than directly calling MPI collectives --- include/ygm/random/alias_table.hpp | 8 ++------ test/test_alias_table.cpp | 6 +++--- 2 files changed, 5 insertions(+), 9 deletions(-) diff --git a/include/ygm/random/alias_table.hpp b/include/ygm/random/alias_table.hpp index c9de59e5..581ffd8e 100644 --- a/include/ygm/random/alias_table.hpp +++ b/include/ygm/random/alias_table.hpp @@ -9,7 +9,6 @@ #include #include #include -#include #include namespace ygm::random { @@ -108,15 +107,12 @@ class alias_table { } void balance_weight() { - MPI_Comm comm = m_comm.get_mpi_comm(); double local_weight = 0.0; for (uint32_t i = 0; i < m_local_items.size(); i++) { local_weight += m_local_items[i].weight; } - double global_weight = 0; - MPI_Allreduce(&local_weight, &global_weight, 1, MPI_DOUBLE, MPI_SUM, comm); - double prfx_sum_weight = 0; - MPI_Exscan(&local_weight, &prfx_sum_weight, 1, MPI_DOUBLE, MPI_SUM, comm); + double global_weight = ygm::sum(local_weight, m_comm); + double prfx_sum_weight = ygm::prefix_sum(local_weight, m_comm); // target_weight = Amount of weight each rank should have after balancing double target_weight = global_weight / m_comm.size(); diff --git a/test/test_alias_table.cpp b/test/test_alias_table.cpp index 1b0b8a51..6bc821b1 100644 --- a/test/test_alias_table.cpp +++ b/test/test_alias_table.cpp @@ -88,7 +88,7 @@ int main(int argc, char** argv) { { // Testing uniform weight distribution std::uniform_int_distribution max_item_weight_dist(50, 100); for (uint32_t i = 0; i < alias_tables_to_construct; i++) { - world.cout0("uniform distribution alias table ", i, " of ", alias_tables_to_construct); + // world.cout0("uniform distribution alias table ", i, " of ", alias_tables_to_construct); ygm::container::map map_of_items(world); uint32_t max_item_weight = max_item_weight_dist(ygm_rng); std::uniform_real_distribution weight_dist(0, max_item_weight); @@ -106,7 +106,7 @@ int main(int argc, char** argv) { std::uniform_int_distribution mean_dist(50, 100); std::uniform_int_distribution std_dev_dist(5, 20); for (uint32_t i = 0; i < alias_tables_to_construct; i++) { - world.cout0("Normal distribution alias table ", i, " of ", alias_tables_to_construct); + // world.cout0("Normal distribution alias table ", i, " of ", alias_tables_to_construct); ygm::container::map map_of_items(world); uint32_t mean = mean_dist(ygm_rng); uint32_t std_dev = std_dev_dist(ygm_rng); @@ -125,7 +125,7 @@ int main(int argc, char** argv) { std::uniform_real_distribution alpha_dist(0.1, 10); std::uniform_real_distribution theta_dist(10, 100); for (uint32_t i = 0; i < alias_tables_to_construct; i++) { - world.cout0("Gamma distribution alias table ", i, " of ", alias_tables_to_construct); + // world.cout0("Gamma distribution alias table ", i, " of ", alias_tables_to_construct); ygm::container::map map_of_items(world); double alpha = alpha_dist(ygm_rng); double theta = theta_dist(ygm_rng);