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

Use Eigen::SelfAdjointView for mults.

Browse Source
This commit is contained in:
David Allemang
2021-10-31 21:19:46 -04:00
parent ec4c1d213c
commit d67768d85c
6 changed files with 93 additions and 59 deletions
Download Patch File Download Diff File Expand all files Collapse all files

3
.gitmodules vendored Normal file
View File

@@ -0,0 +1,3 @@
[submodule "ext/eigen"]
path = ext/eigen
url = https://gitlab.com/libeigen/eigen.git

6
CMakeLists.txt
View File

@@ -3,8 +3,10 @@ project(toddcox-faster)
option(TC_BUILD_EXAMPLE "Build example executables" OFF)
add_library(tc INTERFACE
)
add_subdirectory(ext)
add_library(tc INTERFACE)
target_link_libraries(tc INTERFACE eigen)
target_include_directories(tc INTERFACE include)

2
ext/CMakeLists.txt Normal file
View File

@@ -0,0 +1,2 @@
add_library(eigen INTERFACE)
target_include_directories(eigen INTERFACE eigen)

1
ext/eigen Submodule

Submodule ext/eigen added at b3bea43a2d

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

@@ -10,79 +10,63 @@
#include "rel.hpp"
#include "cosets.hpp"
#include <Eigen/Eigen>
#include <iostream>
namespace tc {
struct Group;
struct SubGroup;
struct Group {
int ngens;
std::vector<std::vector<int>> _mults;
std::string name;
using Matrix = Eigen::MatrixXi;
Group(const Group &) = default;
int ngens;
std::string name;
Matrix _data;
Eigen::SelfAdjointView<Matrix, Eigen::Upper> _mults;
Group(const Group &g)
: ngens(g.ngens),
name(g.name),
_data(g._data),
_mults(_data) {
}
Group(Group &&g) noexcept
: ngens(g.ngens),
name(std::move(g.name)),
_data(std::move(g._data)),
_mults(_data) {
}
explicit Group(
int ngens,
const std::vector<Rel> &rels = {},
std::string name = "G"
) : ngens(ngens), name(std::move(name)) {
_mults.resize(ngens);
for (auto &mult: _mults) {
mult.resize(ngens, 2);
}
for (const auto &rel: rels) {
set(rel);
}
) : ngens(ngens),
name(std::move(name)),
_data(ngens, ngens),
_mults(_data) {
_data.fill(2);
}
void set(const Rel &r) {
_mults[r.gens[0]][r.gens[1]] = r.mult;
_mults[r.gens[1]][r.gens[0]] = r.mult;
Matrix::Scalar &operator()(int a, int b) {
return _mults(a, b);
}
[[nodiscard]] int get(int a, int b) const {
return _mults[a][b];
Matrix::Scalar operator()(int a, int b) const {
return _mults(a, b);
}
[[nodiscard]] std::vector<Rel> get_rels() const {
std::vector<Rel> res;
for (int i = 0; i < ngens - 1; ++i) {
for (int j = i + 1; j < ngens; ++j) {
res.emplace_back(i, j, get(i, j));
res.emplace_back(i, j, operator()(i, j));
}
}
return res;
}
[[nodiscard]] Group product(const Group &other) const {
std::stringstream ss;
ss << name << "*" << other.name;
Group g(ngens + other.ngens, get_rels(), ss.str());
for (const auto &rel: other.get_rels()) {
g.set(rel.shift(ngens));
}
return g;
}
[[nodiscard]] Group power(int p) const {
std::stringstream ss;
ss << name << "^" << p;
Group g(ngens * p, {}, ss.str());
for (const auto &rel: get_rels()) {
for (int off = 0; off < g.ngens; off += ngens) {
g.set(rel.shift(off));
}
}
return g;
}
[[nodiscard]] SubGroup subgroup(
const std::vector<int> &gens
) const;
@@ -103,8 +87,8 @@ namespace tc {
for (size_t i = 0; i < gen_map.size(); ++i) {
for (size_t j = 0; j < gen_map.size(); ++j) {
int mult = parent.get(gen_map[i], gen_map[j]);
set(Rel(i, j, mult));
int mult = parent(gen_map[i], gen_map[j]);
operator()(i, j) = mult;
}
}
}
@@ -114,12 +98,54 @@ namespace tc {
return SubGroup(*this, gens);
}
Group product(const Group &g, const Group &h) {
std::stringstream ss;
ss << g.name << "*" << h.name;
Group res(g.ngens + h.ngens, ss.str());
int off = 0;
for (int i = 0; i < g.ngens; ++i) {
for (int j = i; j < g.ngens; ++j) {
res(i + off, j + off) = g(i, j);
}
}
off += g.ngens;
for (int i = 0; i < h.ngens; ++i) {
for (int j = i; j < h.ngens; ++j) {
res(i + off, j + off) = h(i, j);
}
}
return res;
}
Group power(const Group &g, int p) {
std::stringstream ss;
ss << g.name << "^" << p;
Group res(g.ngens * p, ss.str());
for (int i = 0; i < g.ngens; ++i) {
for (int j = i; j < g.ngens; ++j) {
for (int k = 0; k < p; ++k) {
int off = k * g.ngens;
res(i + off, j + off) = g(i, j);
}
}
}
return res;
}
Group operator*(const Group &g, const Group &h) {
return g.product(h);
return product(g, h);
}
Group operator^(const Group &g, int p) {
return g.power(p);
return power(g, p);
}
}

10
include/tc/groups.hpp
View File

@@ -10,10 +10,10 @@ namespace tc {
Group schlafli(const std::vector<int> &mults, const std::string &name) {
int ngens = (int) mults.size() + 1;
Group g(ngens, {}, name);
Group g(ngens, name);
for (int i = 0; i < (int) mults.size(); i++) {
g.set(Rel(i, i + 1, mults[i]));
g(i, i + 1) = mults[i];
}
return g;
@@ -46,7 +46,7 @@ namespace tc {
ss << "A(" << dim << ")";
if (dim == 0)
return Group(0, {}, ss.str());
return Group(0, ss.str());
const std::vector<int> &mults = std::vector<int>(dim - 1, 3);
@@ -77,7 +77,7 @@ namespace tc {
mults[dim - 2] = 2;
Group g = schlafli(mults, ss.str());
g.set(Rel(1, dim - 1, 3));
g(1, dim - 1) = 3;
return g;
}
@@ -93,7 +93,7 @@ namespace tc {
mults[dim - 2] = 2;
Group g = schlafli(mults, ss.str());
g.set(Rel(2, dim - 1, 3));
g(2, dim - 1) = 3;
return g;
}