me/toddcox-faster
1
0
Fork 0
forked from mirror/toddcox-faster
Code Pull Requests Activity

Templatize tc::Group and tc::Tables

Browse Source
This commit is contained in:
David Allemang
2021-11-01 08:58:11 -04:00
parent 6874faac1d
commit 68065f4a56
6 changed files with 183 additions and 186 deletions
Download Patch File Download Diff File Expand all files Collapse all files

18
example/bench.cpp
View File

@@ -15,24 +15,24 @@ void test(const G &group) {
int order = cosets.size(); int order = cosets.size();
std::cout std::cout
<< std::setw(7) << group.name << ", " << std::setw(7) << group.name() << ", "
<< std::setw(7) << order << ", " << std::setw(7) << order << ", "
<< std::fixed << std::setprecision(6) << diff << "s" << std::fixed << std::setprecision(6) << diff << "s"
<< std::endl; << std::endl;
} }
int main() { int main() {
test(tc::group::H(2)); test(tc::group::H<2>());
test(tc::group::H(3)); test(tc::group::H<3>());
test(tc::group::H(4)); test(tc::group::H<4>());
test(tc::group::T(100)); test(tc::group::T(100));
test(tc::group::T(500)); test(tc::group::T(500));
test(tc::group::T(1000)); test(tc::group::T(1000));
test(tc::group::E(6)); test(tc::group::E<6>());
test(tc::group::E(7)); test(tc::group::E<7>());
test(tc::group::B(6)); test(tc::group::B<6>());
test(tc::group::B(7)); test(tc::group::B<7>());
test(tc::group::B(8)); test(tc::group::B<8>());
return 0; return 0;
} }

1
include/tc/core.hpp
View File

@@ -2,7 +2,6 @@
#include "tc/core/action.hpp" #include "tc/core/action.hpp"
#include "tc/core/path.hpp" #include "tc/core/path.hpp"
#include "tc/core/cosets.hpp"
#include "tc/core/rel.hpp" #include "tc/core/rel.hpp"
#include "tc/core/group.hpp" #include "tc/core/group.hpp"
#include "tc/core/solver.hpp" #include "tc/core/solver.hpp"

39
include/tc/core/cosets.hpp
View File

@@ -1,39 +0,0 @@
#pragma once
#include <vector>
#include "path.hpp"
struct Cosets {
int ngens;
std::vector<int> data;
Path path;
Cosets(const Cosets &) = default;
explicit Cosets(int ngens)
: ngens(ngens) {
}
void add_row() {
data.resize(data.size() + ngens, -1);
path.add_row();
}
void put(int coset, int gen, int target) {
data[coset * ngens + gen] = target;
data[target * ngens + gen] = coset;
if (path.get(target).from_idx == -1) {
path.put(coset, gen, target);
}
}
[[nodiscard]] int get(int coset, int gen) const {
return data[coset * ngens + gen];
}
[[nodiscard]] size_t size() const {
return path.size();
}
};

126
include/tc/core/group.hpp
View File

@@ -8,126 +8,112 @@
#include <queue> #include <queue>
#include "rel.hpp" #include "rel.hpp"
#include "cosets.hpp"
#include <Eigen/Eigen> #include <Eigen/Eigen>
#include <iostream> #include <iostream>
namespace tc { namespace tc {
template<unsigned int Rank>
struct Group; struct Group;
template<unsigned int Rank, unsigned int PRank>
struct SubGroup; struct SubGroup;
struct Group { template<unsigned int Rank>
using Matrix = Eigen::MatrixXi; class Group {
public:
using Matrix = Eigen::Matrix<unsigned int, Rank, Rank>;
int ngens; private:
std::string name; public:
std::string _name;
Matrix _data; Matrix _data;
Eigen::SelfAdjointView<Matrix, Eigen::Upper> _mults; Eigen::SelfAdjointView<Matrix, Eigen::Upper> _mults;
Group(const Group &g) public:
: ngens(g.ngens), Group(const Group<Rank> &g) :
name(g.name), _name(g._name),
_data(g._data), _data(g._data),
_mults(_data) { _mults(_data) {
} }
Group(Group &&g) noexcept Group(Group &&g) noexcept:
: ngens(g.ngens), _name(std::move(g._name)),
name(std::move(g.name)),
_data(std::move(g._data)), _data(std::move(g._data)),
_mults(_data) { _mults(_data) {
} }
explicit Group( explicit Group(std::string name = "G") :
int ngens, _name(std::move(name)),
std::string name = "G" _data(),
) : ngens(ngens),
name(std::move(name)),
_data(ngens, ngens),
_mults(_data) { _mults(_data) {
_data.fill(2); _data.fill(2);
} }
Matrix::Scalar &operator()(int a, int b) { unsigned int rank() const {
return _data.rows();
}
std::string name() const {
return _name;
}
typename Matrix::Scalar &operator()(int a, int b) {
return _mults(a, b); return _mults(a, b);
} }
Matrix::Scalar operator()(int a, int b) const { typename Matrix::Scalar operator()(int a, int b) const {
return _mults(a, b); return _mults(a, b);
} }
[[nodiscard]] std::vector<Rel> get_rels() const { [[nodiscard]] std::vector<Rel> get_rels() const {
std::vector<Rel> res; std::vector<Rel> res;
for (int i = 0; i < ngens - 1; ++i) { for (int i = 0; i < Rank - 1; ++i) {
for (int j = i + 1; j < ngens; ++j) { for (int j = i + 1; j < Rank; ++j) {
res.emplace_back(i, j, operator()(i, j)); res.emplace_back(i, j, _mults(i, j));
} }
} }
return res; return res;
} }
[[nodiscard]] SubGroup subgroup(
const std::vector<int> &gens
) const;
}; };
struct SubGroup : public Group { template<unsigned int GR, unsigned int HR>
std::vector<int> gen_map; Group<GR + HR> product(const Group<GR> &g, const Group<HR> &h) {
const Group &parent;
SubGroup(const Group &parent, std::vector<int> gen_map) : Group(gen_map.size()), parent(parent) {
std::sort(gen_map.begin(), gen_map.end());
this->gen_map = gen_map;
for (size_t i = 0; i < gen_map.size(); ++i) {
for (size_t j = 0; j < gen_map.size(); ++j) {
int mult = parent(gen_map[i], gen_map[j]);
operator()(i, j) = mult;
}
}
}
};
SubGroup Group::subgroup(const std::vector<int> &gens) const {
return SubGroup(*this, gens);
}
Group product(const Group &g, const Group &h) {
std::stringstream ss; std::stringstream ss;
ss << g.name << "*" << h.name; ss << g.name << "*" << h.name;
Group res(g.ngens + h.ngens, ss.str()); Group<GR + HR> res(ss.str());
int off = 0; int off = 0;
for (int i = 0; i < GR; ++i) {
for (int i = 0; i < g.ngens; ++i) { for (int j = i; j < GR; ++j) {
for (int j = i; j < g.ngens; ++j) {
res(i + off, j + off) = g(i, j); res(i + off, j + off) = g(i, j);
} }
} }
off += g.ngens; off += GR;
for (int i = 0; i < h.ngens; ++i) { for (int i = 0; i < HR; ++i) {
for (int j = i; j < h.ngens; ++j) { for (int j = i; j < HR; ++j) {
res(i + off, j + off) = h(i, j); res(i + off, j + off) = h(i, j);
} }
} }
off += HR;
return res; return res;
} }
Group power(const Group &g, int p) { template<unsigned int GR, unsigned int P>
Group<GR * P> power(const Group<GR> &g) {
std::stringstream ss; std::stringstream ss;
ss << g.name << "^" << p; ss << g.name << "^" << P;
Group res(g.ngens * p, ss.str()); Group<GR * P> res(ss.str());
for (int i = 0; i < g.ngens; ++i) { for (int k = 0; k < P; ++k) {
for (int j = i; j < g.ngens; ++j) { int off = k * GR;
for (int k = 0; k < p; ++k) {
int off = k * g.ngens; for (int i = 0; i < GR; ++i) {
for (int j = i; j < GR; ++j) {
res(i + off, j + off) = g(i, j); res(i + off, j + off) = g(i, j);
} }
} }
@@ -136,11 +122,11 @@ namespace tc {
return res; return res;
} }
Group operator*(const Group &g, const Group &h) { // Group operator*(const Group &g, const Group &h) {
return product(g, h); // return product(g, h);
} // }
//
Group operator^(const Group &g, int p) { // Group operator^(const Group &g, int p) {
return power(g, p); // return power(g, p);
} // }
} }

83
include/tc/core/solver.hpp
View File

@@ -3,6 +3,42 @@
#include <memory> #include <memory>
#include <vector> #include <vector>
namespace tc {
struct Cosets {
int ngens;
std::vector<int> data;
Path path;
Cosets(const Cosets &) = default;
explicit Cosets(int ngens)
: ngens(ngens) {
}
void add_row() {
data.resize(data.size() + ngens, -1);
path.add_row();
}
void put(int coset, int gen, int target) {
data[coset * ngens + gen] = target;
data[target * ngens + gen] = coset;
if (path.get(target).from_idx == -1) {
path.put(coset, gen, target);
}
}
[[nodiscard]] int get(int coset, int gen) const {
return data[coset * ngens + gen];
}
[[nodiscard]] size_t size() const {
return path.size();
}
};
}
namespace { namespace {
struct Row { struct Row {
int gnr; int gnr;
@@ -20,20 +56,24 @@ namespace {
} }
}; };
struct Tables { template<unsigned int Rank>
class Tables {
public:
static constexpr unsigned int Rels = Rank * (Rank + 1) / 2 - Rank;
private: private:
std::shared_ptr<int> null_lst_ptr = std::make_shared<int>(); std::shared_ptr<int> null_lst_ptr = std::make_shared<int>();
int ngens;
std::vector<std::shared_ptr<Table>> tables; std::array<std::shared_ptr<Table>, Rels> tables;
std::vector<std::vector<std::shared_ptr<Table>>> deps; std::array<std::vector<std::shared_ptr<Table>>, Rank> deps;
public: public:
explicit Tables(const tc::Group &group) explicit Tables(const tc::Group<Rank> &group) {
: ngens(group.ngens), deps(ngens) { const auto &rels = group.get_rels();
for (const auto &rel: group.get_rels()) { for (int i = 0; i < Rels; ++i) {
const auto &rel = rels[i];
auto table = std::make_shared<Table>(rel); auto table = std::make_shared<Table>(rel);
tables.push_back(table); tables[i] = table;
deps[rel.gens[0]].push_back(table); deps[rel.gens[0]].push_back(table);
deps[rel.gens[1]].push_back(table); deps[rel.gens[1]].push_back(table);
} }
@@ -46,7 +86,7 @@ namespace {
} }
} }
void initialize(int target, const Cosets &cosets) { void initialize(int target, const tc::Cosets &cosets) {
for (auto &table: tables) { for (auto &table: tables) {
const Rel &rel = table->rel; const Rel &rel = table->rel;
Row &row = table->rows[target]; Row &row = table->rows[target];
@@ -64,7 +104,7 @@ namespace {
} }
} }
void learn(int coset, int gen, int target, const Cosets &cosets, std::priority_queue<int> &facts) { void learn(int coset, int gen, int target, const tc::Cosets &cosets, std::priority_queue<int> &facts) {
if (target == coset) { if (target == coset) {
for (auto &table: deps[gen]) { for (auto &table: deps[gen]) {
Row &target_row = table->rows[target]; Row &target_row = table->rows[target];
@@ -92,11 +132,11 @@ namespace {
// forward learn // forward learn
int lst = *target_row.lst; int lst = *target_row.lst;
int gen_ = rel.gens[rel.gens[0] == gen]; int gen_ = rel.gens[rel.gens[0] == gen];
facts.push(lst * ngens + gen_); facts.push(lst * Rank + gen_);
} else if (target_row.gnr == -rel.mult) { } else if (target_row.gnr == -rel.mult) {
// stationary learn // stationary learn
int gen_ = rel.gens[rel.gens[0] == gen]; int gen_ = rel.gens[rel.gens[0] == gen];
facts.push(target * ngens + gen_); facts.push(target * Rank + gen_);
} else if (target_row.gnr == rel.mult - 1) { } else if (target_row.gnr == rel.mult - 1) {
// determined family // determined family
*target_row.lst = target; *target_row.lst = target;
@@ -108,28 +148,29 @@ namespace {
} }
namespace tc { namespace tc {
Cosets solve(const Group &g, const std::vector<int> &sub_gens = {}) { template<unsigned int Rank>
Cosets cosets(g.ngens); tc::Cosets solve(const Group <Rank> &g, const std::vector<int> &sub_gens = {}) {
tc::Cosets cosets(Rank);
cosets.add_row(); cosets.add_row();
if (g.ngens == 0) { if (Rank == 0) {
return cosets; return cosets;
} }
for (int gen: sub_gens) { for (int gen: sub_gens) {
if (gen < g.ngens) if (gen < Rank)
cosets.put(0, gen, 0); cosets.put(0, gen, 0);
} }
Tables tables(g); Tables<Rank> tables(g);
tables.add_row(); tables.add_row();
tables.initialize(0, cosets); tables.initialize(0, cosets);
std::priority_queue<int> facts; std::priority_queue<int> facts;
for (int idx = 0; idx < cosets.data.size(); idx++) { for (int idx = 0; idx < cosets.data.size(); idx++) {
int coset = idx / g.ngens; int coset = idx / Rank;
int gen = idx % g.ngens; int gen = idx % Rank;
if (cosets.get(coset, gen) >= 0) continue; if (cosets.get(coset, gen) >= 0) continue;
int target = cosets.size(); int target = cosets.size();
@@ -147,8 +188,8 @@ namespace tc {
int fact_idx = facts.top(); int fact_idx = facts.top();
facts.pop(); facts.pop();
int coset_ = fact_idx / g.ngens; int coset_ = fact_idx / Rank;
int gen_ = fact_idx % g.ngens; int gen_ = fact_idx % Rank;
if (cosets.get(coset_, gen_) != -1) if (cosets.get(coset_, gen_) != -1)
continue; continue;

102
include/tc/groups.hpp
View File

@@ -7,12 +7,11 @@ namespace tc {
* Construct a group from a (simplified) Schlafli Symbol of the form [a, b, ..., c] * Construct a group from a (simplified) Schlafli Symbol of the form [a, b, ..., c]
* @param mults: The sequence of multiplicites between adjacent generators. * @param mults: The sequence of multiplicites between adjacent generators.
*/ */
Group schlafli(const std::vector<int> &mults, const std::string &name) { template<unsigned int Rank>
int ngens = (int) mults.size() + 1; Group <Rank> schlafli(const std::array<unsigned int, Rank - 1> &mults, const std::string &name) {
Group<Rank> g(name);
Group g(ngens, name); for (int i = 0; i < Rank - 1; i++) {
for (int i = 0; i < (int) mults.size(); i++) {
g(i, i + 1) = mults[i]; g(i, i + 1) = mults[i];
} }
@@ -23,61 +22,68 @@ namespace tc {
* Construct a group from a (simplified) Schlafli Symbol of the form [a, b, ..., c] * Construct a group from a (simplified) Schlafli Symbol of the form [a, b, ..., c]
* @param mults: The sequence of multiplicites between adjacent generators. * @param mults: The sequence of multiplicites between adjacent generators.
*/ */
Group schlafli(const std::vector<int> &mults) { template<unsigned int Rank>
Group <Rank> schlafli(const std::array<unsigned int, Rank - 1> &mults) {
std::stringstream ss; std::stringstream ss;
ss << "["; ss << "[";
if (!mults.empty()) { if (Rank) {
for (size_t i = 0; i < mults.size() - 1; ++i) { for (size_t i = 0; i < Rank - 2; ++i) {
ss << mults[i] << ","; ss << mults[i] << ",";
} }
ss << mults.back(); ss << mults[Rank - 1];
} }
ss << "]"; ss << "]";
return schlafli(mults, ss.str()); return schlafli<Rank>(mults, ss.str());
} }
namespace group { namespace group {
/** /**
* Simplex * Simplex
*/ */
Group A(int dim) { template<unsigned int Rank>
Group <Rank> A() {
std::stringstream ss; std::stringstream ss;
ss << "A(" << dim << ")"; ss << "A(" << Rank << ")";
if (dim == 0) if (Rank == 0)
return Group(0, ss.str()); return Group<Rank>(ss.str());
const std::vector<int> &mults = std::vector<int>(dim - 1, 3); std::array<unsigned int, Rank - 1> mults;
mults.fill(3);
return schlafli(mults, ss.str()); return schlafli<Rank>(mults, ss.str());
} }
/** /**
* Cube, Orthoplex * Cube, Orthoplex
*/ */
Group B(int dim) { template<unsigned int Rank>
Group <Rank> B() {
std::stringstream ss; std::stringstream ss;
ss << "B(" << dim << ")"; ss << "B(" << Rank << ")";
std::vector<int> mults(dim - 1, 3); std::array<unsigned int, Rank - 1> mults;
mults.fill(3);
mults[0] = 4; mults[0] = 4;
return schlafli(mults, ss.str()); return schlafli<Rank>(mults, ss.str());
} }
/** /**
* Demicube, Orthoplex * Demicube, Orthoplex
*/ */
Group D(int dim) { template<unsigned int Rank>
Group <Rank> D() {
std::stringstream ss; std::stringstream ss;
ss << "D(" << dim << ")"; ss << "D(" << Rank << ")";
std::vector<int> mults(dim - 1, 3); std::array<unsigned int, Rank - 1> mults;
mults[dim - 2] = 2; mults.fill(3);
mults[Rank - 2] = 2;
Group g = schlafli(mults, ss.str()); Group<Rank> g = schlafli<Rank>(mults, ss.str());
g(1, dim - 1) = 3; g(1, Rank - 1) = 3;
return g; return g;
} }
@@ -85,15 +91,17 @@ namespace tc {
/** /**
* E groups * E groups
*/ */
Group E(int dim) { template<unsigned int Rank>
Group <Rank> E() {
std::stringstream ss; std::stringstream ss;
ss << "E(" << dim << ")"; ss << "E(" << Rank << ")";
std::vector<int> mults(dim - 1, 3); std::array<unsigned int, Rank - 1> mults;
mults[dim - 2] = 2; mults.fill(3);
mults[Rank - 2] = 2;
Group g = schlafli(mults, ss.str()); Group<Rank> g = schlafli<Rank>(mults, ss.str());
g(2, dim - 1) = 3; g(2, Rank - 1) = 3;
return g; return g;
} }
@@ -101,58 +109,60 @@ namespace tc {
/** /**
* 24 Cell * 24 Cell
*/ */
Group F4() { Group<4> F4() {
return schlafli({3, 4, 3}, "F4"); return schlafli<4>({3, 4, 3}, "F4");
} }
/** /**
* Hexagon * Hexagon
*/ */
Group G2() { Group<2> G2() {
return schlafli({6}, "G2"); return schlafli<2>({6}, "G2");
} }
/** /**
* Icosahedron * Icosahedron
*/ */
Group H(int dim) { template<unsigned int Rank>
Group <Rank> H() {
std::stringstream ss; std::stringstream ss;
ss << "H(" << dim << ")"; ss << "H(" << Rank << ")";
std::vector<int> mults(dim - 1, 3); std::array<unsigned int, Rank - 1> mults;
mults.fill(3);
mults[0] = 5; mults[0] = 5;
return schlafli(mults, ss.str()); return schlafli<Rank>(mults, ss.str());
} }
/** /**
* Polygonal * Polygonal
*/ */
Group I2(int n) { Group<2> I2(unsigned int n) {
std::stringstream ss; std::stringstream ss;
ss << "I2(" << n << ")"; ss << "I2(" << n << ")";
return schlafli({n}, ss.str()); return schlafli<2>({n}, ss.str());
} }
/** /**
* Toroidal. I2(n) * I2(m) * Toroidal. I2(n) * I2(m)
*/ */
Group T(int n, int m) { Group<4> T(unsigned int n, unsigned int m) {
std::stringstream ss; std::stringstream ss;
ss << "T(" << n << "," << m << ")"; ss << "T(" << n << "," << m << ")";
return schlafli({n, 2, m}, ss.str()); return schlafli<4>({n, 2, m}, ss.str());
} }
/** /**
* Toroidal. T(n, n) * Toroidal. T(n, n)
*/ */
Group T(int n) { Group<4> T(unsigned int n) {
std::stringstream ss; std::stringstream ss;
ss << "T(" << n << ")"; ss << "T(" << n << ")";
return schlafli({n, 2, n}, ss.str()); return schlafli<4>({n, 2, n}, ss.str());
} }
} }
} }