From a1b0ae31dc6e427227598da66afc5a10293091b2 Mon Sep 17 00:00:00 2001 From: David Allemang Date: Fri, 12 Nov 2021 22:14:46 -0500 Subject: [PATCH] Un-template complex solver --- examples/complex.cpp | 27 ++++--- include/tc/complex.hpp | 172 ++++++++++++++++------------------------- 2 files changed, 78 insertions(+), 121 deletions(-) diff --git a/examples/complex.cpp b/examples/complex.cpp index b061a81..59b4780 100644 --- a/examples/complex.cpp +++ b/examples/complex.cpp @@ -5,11 +5,11 @@ #include #include -template -void test(const G &group) { +template +void test(const G &group, unsigned int N) { auto s = std::clock(); auto combos = tc::combinations(group.gens, N - 1); - auto data = tc::merge(tc::hull(group, combos, {})); + auto data = tc::merge(tc::hull(group, combos, {})); auto e = std::clock(); double diff = (double) (e - s) / CLOCKS_PER_SEC; @@ -24,17 +24,16 @@ void test(const G &group) { } int main() { - test<4>(tc::group::H(4)); - test<4>(tc::group::B(4)); - test<4>(tc::group::B(5)); - test<4>(tc::group::B(6)); - test<4>(tc::group::E(6)); - - test<3>(tc::group::H(3)); - test<3>(tc::group::B(4)); - test<3>(tc::group::B(5)); - test<3>(tc::group::B(6)); - test<3>(tc::group::E(6)); + test(tc::group::H(4), 4); + test(tc::group::B(4), 4); + test(tc::group::B(5), 4); + test(tc::group::B(6), 4); + test(tc::group::E(6), 4); + test(tc::group::H(3), 3); + test(tc::group::B(4), 3); + test(tc::group::B(5), 3); + test(tc::group::B(6), 3); + test(tc::group::E(6), 3); return 0; } diff --git a/include/tc/complex.hpp b/include/tc/complex.hpp index cf6028f..240f7f7 100644 --- a/include/tc/complex.hpp +++ b/include/tc/complex.hpp @@ -6,6 +6,7 @@ #include #include #include +#include namespace tc { std::vector combinations(const Symbol &symbol, size_t srank) { @@ -30,38 +31,41 @@ namespace tc { return combos; } -/** - * An primitive stage N indices. - * @tparam N - */ - template - struct Primitive { - static_assert(N > 0, "Primitives must contain at least one point. Primitive<0> or lower is impossible."); + Symbol fan(const Symbol &prim, unsigned root) { + Symbol res(prim.size() + 1); + res << prim, root; + return res; + } - std::array inds; - - Primitive() = default; - - Primitive(const Primitive &) = default; - - Primitive(const Primitive &sub, unsigned root) { - std::copy(sub.inds.begin(), sub.inds.end(), inds.begin()); - inds[N - 1] = root; + std::vector fan(const std::vector& prims, int root) { + std::vector res; + res.reserve(prims.size()); + for (const auto &prim: prims) { + Symbol s(prim.size() + 1); + s << prim, root; + res.push_back(s); } + return res; + } - ~Primitive() = default; - - inline void flip() { - if (N > 1) std::swap(inds[0], inds[1]); + void flip(Symbol &prim) { + if (prim.size() > 1) { + std::swap(prim[0], prim[1]); } + } - void apply(const tc::Cosets &table, unsigned int gen) { - for (auto &ind: inds) { - ind = table.get(ind, gen); - } - flip(); + void apply(const tc::Cosets &table, unsigned int gen, Symbol &prim) { + for (auto &ind: prim) { + ind = table.get(ind, gen); } - }; + flip(prim); + } + + void apply(const tc::Cosets &table, unsigned int gen, std::vector &prims) { + for (auto &prim: prims) { + apply(table, gen, prim); + } + } /** * Produce a list of all generators for the group context. The range [0..group.rank). @@ -94,35 +98,20 @@ namespace tc { return i & 1; } - -/** - * Apply some context transformation to all primitives of this mesh. - */ - template - std::vector> apply(std::vector> prims, const tc::Cosets &table, unsigned int gen) { - for (auto &prim: prims) { - prim.apply(table, gen); - } - return prims; - } - /** * Reverse the orientation of all primitives in this mesh. */ - template - void flip(std::vector> prims) { + void flip(std::vector &prims) { for (auto &prim: prims) { - prim.flip(); + flip(prim); } } /** * Convert the indexes of this mesh to those of a different context, using g_gens to build the parent context and sg_gens to build this context. */ - template - [[nodiscard]] - std::vector> recontext( - std::vector> prims, + void recontext( + std::vector &prims, const tc::Group &context, const Symbol &g_gens, const Symbol &sg_gens @@ -136,30 +125,26 @@ namespace tc { return table.get(coset, gen); }); - std::vector> res(prims); - for (Primitive &prim: res) { - for (auto &ind: prim.inds) { + for (Symbol &prim: prims) { + for (auto &ind: prim) { ind = map[ind]; } } if (get_parity(context, g_gens, sg_gens) == 1) - flip(res); - - return res; + flip(prims); } /** * Union several meshes of the same dimension */ - template - std::vector> merge(const std::vector>> &meshes) { + std::vector merge(const std::vector> &meshes) { size_t size = 0; for (const auto &mesh: meshes) { size += mesh.size(); } - std::vector> res; + std::vector res; res.reserve(size); for (const auto &mesh: meshes) { res.insert(res.end(), mesh.begin(), mesh.end()); @@ -168,100 +153,73 @@ namespace tc { return res; } - template - [[nodiscard]] - std::vector>> each_tile( - std::vector> prims, + std::vector> each_tile( + std::vector base, const tc::Group &context, const Symbol &g_gens, const Symbol &sg_gens ) { - std::vector> base = recontext(prims, context, g_gens, sg_gens); + recontext(base, context, g_gens, sg_gens); const auto table = solve(context, g_gens, Symbol(0)); const auto path = solve(context, g_gens, sg_gens).path(); auto _gens = generators(context); - auto res = path.walk(base, _gens, [&table](auto from, auto gen){ - return apply(from, table, gen); + auto res = path.walk(base, _gens, [&table](auto from, auto &gen) { + apply(table, gen, from); + return from; }); return res; } - template [[nodiscard]] - std::vector> tile( - std::vector> prims, + std::vector tile( + std::vector base, const tc::Group &context, const Symbol &g_gens, const Symbol &sg_gens ) { - auto res = each_tile(prims, context, g_gens, sg_gens); + auto res = each_tile(std::move(base), context, g_gens, sg_gens); return merge(res); } -/** - * Produce a mesh of higher dimension by fanning a single point to all primitives in this mesh. - */ - template - [[nodiscard]] - std::vector> fan(std::vector> prims, int root) { - std::vector> res(prims.size()); - std::transform(prims.begin(), prims.end(), res.begin(), - [root](const Primitive &prim) { - return Primitive(prim, root); - } - ); - return res; - } - /** * Produce a mesh of primitives that fill out the volume of the subgroup generated by generators g_gens within the group context */ - template - std::vector> triangulate( + std::vector triangulate( const tc::Group &context, const Symbol &g_gens ) { - if (g_gens.size() + 1 != N) // todo make static assert - throw std::logic_error("g_gens size must be one less than N"); + if (g_gens.size() == 0) { + std::vector res; + Symbol prims(1); + prims.setZero(); + res.push_back(prims); + return res; + } const auto &combos = combinations(g_gens, g_gens.size() - 1); - std::vector>> meshes; + std::vector> meshes; for (const auto &sg_gens: combos) { - auto base = triangulate(context, sg_gens); + auto base = triangulate(context, sg_gens); auto raised = tile(base, context, g_gens, sg_gens); raised.erase(raised.begin(), raised.begin() + base.size()); - meshes.push_back(fan(raised, 0)); + auto fanned = fan(raised, 0); + meshes.push_back(fanned); } - return merge(meshes); + const std::vector &result = merge(meshes); + return result; } -/** - * Single-index primitives should not be further triangulated. - */ - template<> - std::vector> triangulate( - const tc::Group &context, - const Symbol &g_gens - ) { - if (g_gens.size() != 0) // todo make static assert - throw std::logic_error("g_gens must be empty for a trivial Mesh"); - - std::vector> res; - res.emplace_back(); - return res; - } - - template + template auto hull(const tc::Group &group, T all_sg_gens, const std::vector &exclude) { - std::vector>> parts; + std::vector> parts; auto g_gens = group.gens; for (const Symbol &sg_gens: all_sg_gens) { bool excluded = false; @@ -273,7 +231,7 @@ namespace tc { } if (excluded) continue; - const auto &base = triangulate(group, sg_gens); + const auto &base = triangulate(group, sg_gens); const auto &tiles = each_tile(base, group, g_gens, sg_gens); for (const auto &tile: tiles) { parts.push_back(tile);