diff --git a/Testing/common.hpp b/Testing/common.hpp index 165eb41..c4aa908 100644 --- a/Testing/common.hpp +++ b/Testing/common.hpp @@ -21,11 +21,10 @@ std::vector parse_vec(const std::string &part) { return res; } -template -int compute( - const tc::Group &group, +size_t compute( + const tc::Group &group, const std::vector &vgens ) { - auto table = tc::solve(group, vgens); + auto table = tc::solve(group, vgens); return table.order(); } diff --git a/Testing/schlafli.cpp b/Testing/schlafli.cpp index dd990cd..9631970 100644 --- a/Testing/schlafli.cpp +++ b/Testing/schlafli.cpp @@ -7,34 +7,10 @@ int main(int argc, char **argv) { auto vgens = parse_vec(argv[2]); auto target = std::stoul(argv[3]); - int order; - switch (vsymbol.size()) { - case 0: { - Eigen::Vector symbol(vsymbol.data()); - order = compute<1>(tc::schlafli<1>(symbol), vgens); - break; - } - case 1: { - Eigen::Vector symbol(vsymbol.data()); - order = compute<2>(tc::schlafli<2>(symbol), vgens); - break; - } - case 2: { - Eigen::Vector symbol(vsymbol.data()); - order = compute<3>(tc::schlafli<3>(symbol), vgens); - break; - } - case 3: { - Eigen::Vector symbol(vsymbol.data()); - order = compute<4>(tc::schlafli<4>(symbol), vgens); - break; - } - case 4: { - Eigen::Vector symbol(vsymbol.data()); - order = compute<5>(tc::schlafli<5>(symbol), vgens); - break; - } - } + tc::Symbol symbol(vsymbol.size()); + symbol << Eigen::Map(vsymbol.data(), vsymbol.size()); + tc::Group group = tc::schlafli(symbol); + auto order = compute(group, vgens); std::cout << "Order: " << order << ":" << target << std::endl; return order != target; diff --git a/Testing/special.cpp b/Testing/special.cpp index 2dba89e..db62c6f 100644 --- a/Testing/special.cpp +++ b/Testing/special.cpp @@ -9,31 +9,29 @@ int main(int argc, char **argv) { unsigned int order; + tc::Group group(0); + if (name == "E6") { - auto group = tc::group::E<6>(); - order = compute<6>(group, vgens); + group = tc::group::E(6); } if (name == "E7") { - auto group = tc::group::E<7>(); - order = compute<7>(group, vgens); + group = tc::group::E(7); } if (name == "E8") { - auto group = tc::group::E<8>(); - order = compute<8>(group, vgens); + group = tc::group::E(8); } if (name == "B6") { - auto group = tc::group::B<6>(); - order = compute<6>(group, vgens); + group = tc::group::B(6); } if (name == "B7") { - auto group = tc::group::B<7>(); - order = compute<7>(group, vgens); + group = tc::group::B(7); } if (name == "B8") { - auto group = tc::group::B<8>(); - order = compute<8>(group, vgens); + group = tc::group::B(8); } + order = compute(group, vgens); + std::cout << "Order: " << order << ":" << target << std::endl; return order != target; } diff --git a/examples/bench.cpp b/examples/bench.cpp index ee952a2..08904eb 100644 --- a/examples/bench.cpp +++ b/examples/bench.cpp @@ -22,17 +22,17 @@ void test(const G &group) { } int main() { - test(tc::group::H<2>()); - test(tc::group::H<3>()); - test(tc::group::H<4>()); + test(tc::group::H(2)); + test(tc::group::H(3)); + test(tc::group::H(4)); test(tc::group::T(100)); test(tc::group::T(500)); test(tc::group::T(1000)); - test(tc::group::E<6>()); - test(tc::group::E<7>()); - test(tc::group::B<6>()); - test(tc::group::B<7>()); - test(tc::group::B<8>()); + test(tc::group::E(6)); + test(tc::group::E(7)); + test(tc::group::B(6)); + test(tc::group::B(7)); + test(tc::group::B(8)); return 0; } diff --git a/examples/group.cpp b/examples/group.cpp index 8adb148..e197b43 100644 --- a/examples/group.cpp +++ b/examples/group.cpp @@ -4,20 +4,21 @@ #include int main() { - constexpr unsigned int Rank = 5; - constexpr unsigned int SRank = 3; + tc::Symbol symbol(4); + symbol << 5, 3, 2, 3; + tc::Group group = tc::schlafli(symbol); +// tc::Group group = tc::group::E(7); - tc::Group group = tc::schlafli({5, 3, 2, 3}); -// tc::Group group = tc::group::E(); + size_t srank = 3; - tc::SubGroups subs = tc::subgroups(group); + tc::SubGroups subs = tc::subgroups(group, srank); std::cout << "Group " << group.name << " (" << subs.size() << " subgroups)" << std::endl; std::cout << group << std::endl; for (const auto &sub: subs) { - for (int i = 0; i < SRank; ++i) { - for (int j = 0; j < SRank; ++j) { + for (int i = 0; i < srank; ++i) { + for (int j = 0; j < srank; ++j) { auto sub_mult = sub(i, j); auto src_mult = group(sub.gens(i), sub.gens(j)); diff --git a/examples/path.cpp b/examples/path.cpp index f3aca5f..da7c584 100644 --- a/examples/path.cpp +++ b/examples/path.cpp @@ -6,7 +6,7 @@ #include int main() { - auto cube = tc::group::B<3>(); + auto cube = tc::group::B(3); auto vars = tc::solve(cube); std::string start; diff --git a/include/tc/cosets.hpp b/include/tc/cosets.hpp index 4f50b32..fc7002e 100644 --- a/include/tc/cosets.hpp +++ b/include/tc/cosets.hpp @@ -3,58 +3,65 @@ #include namespace tc { - template class Path; - template class Cosets { private: std::vector data; + size_t _rank; public: - Cosets(const Cosets &) = default; + Cosets(const Cosets &) = default; - Cosets() = default; + explicit Cosets(size_t rank) : _rank(rank) {} void add_row() { - data.resize(data.size() + Rank, -1); + data.resize(data.size() + rank(), -1); } void put(int coset, int gen, int target) { - data[coset * Rank + gen] = target; - data[target * Rank + gen] = coset; + data[coset * rank() + gen] = target; + data[target * rank() + gen] = coset; } [[nodiscard]] int get(int coset, int gen) const { - return data[coset * Rank + gen]; + return data[coset * rank() + gen]; + } + + [[nodiscard]] size_t rank() const { + return _rank; } [[nodiscard]] size_t order() const { - return data.size() / Rank; + return data.size() / _rank; } - Path path() const; + Path path() const; }; - template class Path { private: - friend class Cosets; + friend class Cosets; std::vector source; std::vector gen; + size_t _rank; size_t _order; - explicit Path(size_t order) : _order(order), source(order), gen(order) {} + explicit Path(size_t rank, size_t order) : _rank(rank), _order(order), source(order), gen(order) {} public: - size_t order() const { + [[nodiscard]] size_t rank() const { + return _rank; + } + + [[nodiscard]] size_t order() const { return _order; } template - std::vector walk(const T &start, const F &op) { - std::vector res; + std::vector walk(const T &start, const F &op) { + std::vector res; res.reserve(order()); res.push_back(start); @@ -67,20 +74,19 @@ namespace tc { } template - std::vector walk(const T &start, const E &gens, const F &op) { + std::vector walk(const T &start, const E &gens, const F &op) { return walk(start, [&](const T &s, const int g) { return op(s, gens[g]); }); } }; - template - Path Cosets::path() const { - Path res(order()); + Path Cosets::path() const { + Path res(rank(), order()); std::vector set(order()); for (int coset = 0; coset < order(); ++coset) { - for (int gen = 0; gen < Rank; ++gen) { + for (int gen = 0; gen < rank(); ++gen) { int target = get(coset, gen); if (!set[target]) { diff --git a/include/tc/group.hpp b/include/tc/group.hpp index 5a6b17b..7ac586f 100644 --- a/include/tc/group.hpp +++ b/include/tc/group.hpp @@ -16,38 +16,39 @@ namespace { } namespace tc { - template - using Symbol = Eigen::Vector; + using Symbol = Eigen::Vector; - template - Symbol iota() { - Symbol res; - for (int i = 0; i < Rank; ++i) { - res(i) = i; - } - return res; - } + using MatrixXui = Eigen::Matrix; /// A Coxeter Matrix - template - class Group : public Eigen::Matrix { + class Group : public MatrixXui { public: - using Base = Eigen::Matrix; + using Base = MatrixXui; std::string name = "G"; - Symbol gens = iota(); + Symbol gens; - using Base::Base; + explicit Group(size_t rank) : Base(rank, rank), gens(rank) { + for (Eigen::Index i = 0; i < rank; ++i) { + gens(i) = i; + } + } + + [[nodiscard]] size_t rank() const { + return rows(); + } }; - template - Group subgroup(const Group &group, Symbol gens) { - Group res; + Group subgroup(const Group &group, const Symbol &gens) { + size_t rank = group.size(); + size_t srank = gens.size(); + + Group res(srank); res.name = group.name + ":" + stringify(gens); res.gens = gens; - for (int i = 0; i < SRank; ++i) { - for (int j = 0; j < SRank; ++j) { + for (Eigen::Index i = 0; i < srank; ++i) { + for (Eigen::Index j = 0; j < srank; ++j) { res(i, j) = group(gens[i], gens[j]); } } @@ -55,17 +56,17 @@ namespace tc { return res; } - template - Symbol inverse(Symbol gens) { - Symbol res; + Symbol inverse(size_t rank, const Symbol &gens) { + size_t srank = gens.size(); + Symbol res(rank); res.fill(0); - for (int i = 0; i < SRank; ++i) { + for (int i = 0; i < srank; ++i) { res(gens(i)) = i; } return res; } - constexpr unsigned int Factorial(unsigned int n) { + unsigned int factorial(unsigned int n) { unsigned int res = 1; for (int i = 1; i <= n; ++i) { res *= i; @@ -73,28 +74,28 @@ namespace tc { return res; } - constexpr unsigned int Choose(unsigned int n, unsigned int k) { - return Factorial(n) / Factorial(k) / Factorial(n - k); + unsigned int choose(unsigned int n, unsigned int k) { + return factorial(n) / factorial(k) / factorial(n - k); } - template - using SubGroups = std::array, Choose(Rank, SRank)>; + using SubGroups = std::vector; - template - SubGroups subgroups(const Group &group) { - std::vector mask(Rank, false); - std::fill_n(mask.begin(), SRank, true); + SubGroups subgroups(const Group &group, size_t srank) { + size_t rank = group.rank(); - SubGroups res; + std::vector mask(rank, false); + std::fill_n(mask.begin(), srank, true); - size_t i = 0; - Symbol row; + SubGroups res; + res.reserve(choose(rank, srank)); + + Symbol row(srank); do { - for (int j = 0, k = 0; j < Rank; ++j) { + for (int j = 0, k = 0; j < rank; ++j) { if (mask[j]) row(k++) = j; } - res[i++] = subgroup(group, row); + res.push_back(subgroup(group, row)); } while (std::prev_permutation(mask.begin(), mask.end())); return res; @@ -103,15 +104,16 @@ namespace tc { /** * Create a named coxeter matrix from a simplified schlafli symbol */ - template - Group schlafli(const Symbol &mults, const std::string &name) { - Group res; + Group schlafli(const Symbol &mults, const std::string &name) { + size_t rank = mults.size() + 1; + + Group res(rank); res.name = name; res.fill(2); res.diagonal().fill(1); - res.topRightCorner(Rank - 1, Rank - 1).diagonal() << mults; - res.bottomLeftCorner(Rank - 1, Rank - 1).diagonal() << mults; + res.topRightCorner(rank - 1, rank - 1).diagonal() << mults; + res.bottomLeftCorner(rank - 1, rank - 1).diagonal() << mults; return res; } @@ -119,38 +121,36 @@ namespace tc { /** * Create a coxeter matrix from a simplified schlafli symbol. */ - template - Group schlafli(const Symbol &mults) { - return schlafli(mults, stringify(mults)); + Group schlafli(const Symbol &mults) { + return schlafli(mults, stringify(mults)); } - template - Group product(const Group &g, const Group
&h) { - Group res; + Group product(const Group &g, const Group &h) { + Group res(g.rank() + h.rank()); res.name = g.name + "*" + h.name; res.fill(2); - int off = 0; - res.block(off, off, GR, GR) << g.array() + off; - off += GR; + Eigen::Index off = 0; + res.block(off, off, g.rank(), g.rank()) << g.array() + off; + off += (Eigen::Index) g.rank(); - res.block(off, off, HR, HR) << h.array() + off; - off += HR; + res.block(off, off, h.rank(), h.rank()) << h.array() + off; + off += (Eigen::Index) h.rank(); return res; } - template - Group power(const Group &g) { - Group res; - res.name = g.name + "^" + P; + Group power(const Group &g, size_t p) { + Group res(g.rank() * p); + + res.name = g.name + "^" + std::to_string(p); res.fill(2); - for (int k = 0; k < P; ++k) { - int off = k * GR; - res.block(off, off, GR, GR) << g.array() + off; + for (Eigen::Index k = 0; k < p; ++k) { + auto off = (Eigen::Index) g.rank() * k; + res.block(off, off, g.rank(), g.rank()) << g.array() + off; } return res; diff --git a/include/tc/groups.hpp b/include/tc/groups.hpp index 81a2c54..56fe385 100644 --- a/include/tc/groups.hpp +++ b/include/tc/groups.hpp @@ -6,13 +6,11 @@ namespace tc::group { /** * Universal Coxeter Group */ - template - Group U() { - std::stringstream ss; - ss << "U(" << Rank << ")"; + Group U(size_t rank) { + std::string name = "U(" + std::to_string(rank) + ")"; - Group res; - res.name = ss.str(); + Group res(rank); + res.name = name; res.fill(2); return res; } @@ -20,53 +18,47 @@ namespace tc::group { /** * Simplex */ - template - Group A() { - std::stringstream ss; - ss << "A(" << Rank << ")"; + Group A(size_t rank) { + std::string name = "A(" + std::to_string(rank) + ")"; - if (Rank == 0) { - Group res; - res.name = ss.str(); + if (rank == 0) { + Group res(rank); + res.name = name; return res; } - tc::Symbol mults; - mults.fill(3); + tc::Symbol symbol(rank - 1); + symbol.fill(3); - return schlafli(mults, ss.str()); + return schlafli(symbol, name); } /** * Cube, Orthoplex */ - template - Group B() { - std::stringstream ss; - ss << "B(" << Rank << ")"; + Group B(size_t rank) { + std::string name = "B(" + std::to_string(rank) + ")"; - tc::Symbol mults; - mults.fill(3); - mults(0) = 4; + tc::Symbol symbol(rank - 1); + symbol.fill(3); + symbol(0) = 4; - return schlafli(mults, ss.str()); + return schlafli(symbol, name); } /** * Demicube, Orthoplex */ - template - Group D() { - std::stringstream ss; - ss << "D(" << Rank << ")"; + Group D(size_t rank) { + std::string name = "D(" + std::to_string(rank) + ")"; - tc::Symbol mults; - mults.fill(3); - mults(Rank - 2) = 2; + tc::Symbol symbol(rank - 1); + symbol.fill(3); + symbol((Eigen::Index) rank - 2) = 2; - Group g = schlafli(mults, ss.str()); - g(1, Rank - 1) = 3; - g(Rank - 1, 1) = 3; + Group g = schlafli(symbol, name); + g(1, (Eigen::Index) rank - 1) = 3; + g((Eigen::Index) rank - 1, 1) = 3; return g; } @@ -74,18 +66,16 @@ namespace tc::group { /** * E groups */ - template - Group E() { - std::stringstream ss; - ss << "E(" << Rank << ")"; + Group E(size_t rank) { + std::string name = "E(" + std::to_string(rank) + ")"; - tc::Symbol mults; - mults.fill(3); - mults(Rank - 2) = 2; + tc::Symbol symbol(rank - 1); + symbol.fill(3); + symbol((Eigen::Index) rank - 2) = 2; - Group g = schlafli(mults, ss.str()); - g(2, Rank - 1) = 3; - g(Rank - 1, 2) = 3; + Group g = schlafli(symbol, name); + g(2, (Eigen::Index) rank - 1) = 3; + g((Eigen::Index) rank - 1, 2) = 3; return g; } @@ -93,59 +83,66 @@ namespace tc::group { /** * 24 Cell */ - Group<4> F4() { - return schlafli<4>({3, 4, 3}, "F4"); + Group F4() { + tc::Symbol symbol(3); + symbol << 3, 4, 3; + return schlafli(symbol, "F4"); } /** * Hexagon */ - Group<2> G2() { - return schlafli<2>(tc::Symbol<1>(6), "G2"); + Group G2() { + tc::Symbol symbol(1); + symbol << 6; + return schlafli(symbol, "G2"); } /** * Icosahedron */ - template - Group H() { - std::stringstream ss; - ss << "H(" << Rank << ")"; + Group H(size_t rank) { + std::string name = "H(" + std::to_string(rank) + ")"; - tc::Symbol mults; - mults.fill(3); - mults(0) = 5; + tc::Symbol symbol(rank - 1); + symbol.fill(3); + symbol(0) = 5; - return schlafli(mults, ss.str()); + return schlafli(symbol, name); } /** * Polygonal */ - Group<2> I2(unsigned int n) { - std::stringstream ss; - ss << "I2(" << n << ")"; + Group I2(unsigned int n) { + std::string name = "I2(" + std::to_string(n) + ")"; - return schlafli<2>(tc::Symbol<1>(n), ss.str()); + tc::Symbol symbol(1); + symbol << n; + return schlafli(symbol, name); } /** * Toroidal. I2(n) * I2(m) */ - Group<4> T(unsigned int n, unsigned int m) { - std::stringstream ss; - ss << "T(" << n << "," << m << ")"; + Group T(unsigned int n, unsigned int m) { + std::string name = "T(" + std::to_string(n) + "," + std::to_string(m) + ")"; - return schlafli<4>({n, 2, m}, ss.str()); + tc::Symbol symbol(3); + symbol << n, 2, m; + + return schlafli(symbol, name); } /** * Toroidal. T(n, n) */ - Group<4> T(unsigned int n) { - std::stringstream ss; - ss << "T(" << n << ")"; + Group T(unsigned int n) { + std::string name = "T(" + std::to_string(n) + ")"; - return schlafli<4>({n, 2, n}, ss.str()); + tc::Symbol symbol(3); + symbol << n, 2, n; + + return schlafli(symbol, name); } } diff --git a/include/tc/solver.hpp b/include/tc/solver.hpp index a20d6ca..e0bd4b0 100644 --- a/include/tc/solver.hpp +++ b/include/tc/solver.hpp @@ -20,7 +20,7 @@ namespace { public: int i, j, mult; - std::vector rows; + std::vector rows; public: explicit Table(int i, int j, int mult) : @@ -60,31 +60,38 @@ namespace { } }; - template class Tables { - public: - static constexpr unsigned int Rels = Rank * (Rank + 1) / 2 - Rank; - private: int *null_lst_ptr = new int; BlockAllocator alloc; - std::array , Rels> tables; - std::array >, Rank> - deps; + std::vector> tables; + std::vector>> deps; + + size_t _rank; + size_t _rels; public: - explicit Tables(const tc::Group &group) { - for (int i = 0, irel = 0; i < Rank - 1; ++i) { - for (int j = i + 1; j < Rank; ++j, ++irel) { + explicit Tables(const tc::Group &group) : _rank(group.rank()), _rels(rank() * (rank() + 1) / 2 - rank()) { + deps.resize(rank()); + for (int i = 0; i < rank() - 1; ++i) { + for (int j = i + 1; j < rank(); ++j) { auto table = std::make_shared(i, j, group(i, j)); - tables[irel] = table; + tables.push_back(table); deps[i].push_back(table); deps[j].push_back(table); } } } + [[nodiscard]] size_t rank() const { + return _rank; + } + + [[nodiscard]] size_t rels() const { + return _rels; + } + void add_row() { // std::vector already does block allocation. for (const auto &table: tables) { @@ -92,7 +99,7 @@ namespace { } } - void initialize(int target, const tc::Cosets &cosets) { + void initialize(int target, const tc::Cosets &cosets) { for (auto &table: tables) { Row &row = table->rows[target]; @@ -113,7 +120,7 @@ namespace { delete null_lst_ptr; } - void learn(int coset, int gen, int target, const tc::Cosets &cosets, std::priority_queue &facts) { + void learn(int coset, int gen, int target, const tc::Cosets &cosets, std::priority_queue &facts) { if (target == coset) { for (auto &table: deps[gen]) { Row &target_row = table->rows[target]; @@ -140,11 +147,11 @@ namespace { // forward learn int lst = *target_row.lst; int gen_ = (table->i == gen) ? table->j : table->i; - facts.push(lst * Rank + gen_); + facts.push(lst * rank() + gen_); } else if (target_row.gnr == -table->mult) { // stationary learn int gen_ = (table->i == gen) ? table->j : table->i; - facts.push(target * Rank + gen_); + facts.push(target * rank() + gen_); } else if (target_row.gnr == table->mult - 1) { // determined family *target_row.lst = target; @@ -159,35 +166,36 @@ namespace tc { /** * Assumes that g is a coxeter group - that is, self-adjoint and the diagonal is 2. */ - template - tc::Cosets solve(const Group &group, const std::vector &sub_gens = {}) { - tc::Cosets cosets; + tc::Cosets solve(const Group &group, const std::vector &sub_gens = {}) { + size_t rank = group.rank(); + + tc::Cosets cosets(rank); cosets.add_row(); - if (Rank == 0) { + if (rank == 0) { return cosets; } for (unsigned int gen: sub_gens) { - if (gen < Rank) + if (gen < rank) cosets.put(0, gen, 0); } - Tables tables(group); + Tables tables(group); tables.add_row(); tables.initialize(0, cosets); - std::priority_queue facts; + std::priority_queue facts; for (int coset = 0; coset < cosets.order(); coset++) { - for (int gen = 0; gen < Rank; ++gen) { + for (int gen = 0; gen < rank; ++gen) { if (cosets.get(coset, gen) >= 0) continue; // todo vector set int target = cosets.order(); cosets.add_row(); tables.add_row(); - facts.push(coset * Rank + gen); + facts.push(coset * rank + gen); // todo nothing before the current coset will be used. // delete all table rows using old cosets to free memory early. @@ -198,8 +206,8 @@ namespace tc { int fact_idx = facts.top(); facts.pop(); - int coset_ = fact_idx / Rank; - int gen_ = fact_idx % Rank; + int coset_ = fact_idx / rank; + int gen_ = fact_idx % rank; if (cosets.get(coset_, gen_) != -1) continue;