WIP (broken): Refactor components

- move group, root into Structure<>
- semantic type aliases like VectorDf, MatrixGXui, etc.
- separate components and systems files

.gitignore

@@ -1,4 +1,288 @@
-.idea
+### C++ template
+# Prerequisites
+*.d
+
+# Compiled Object files
+*.slo
+*.lo
+*.o
+*.obj
+
+# Precompiled Headers
+*.gch
+*.pch
+
+# Compiled Dynamic libraries
+*.so
+*.dylib
+*.dll
+
+# Fortran module files
+*.mod
+*.smod
+
+# Compiled Static libraries
+*.lai
+*.la
+*.a
+*.lib
+
+# Executables
+*.exe
+*.out
+*.app
+
+### JetBrains template
+# Covers JetBrains IDEs: IntelliJ, RubyMine, PhpStorm, AppCode, PyCharm, CLion, Android Studio, WebStorm and Rider
+# Reference: https://intellij-support.jetbrains.com/hc/en-us/articles/206544839
+
+# User-specific stuff
+.idea/**/workspace.xml
+.idea/**/tasks.xml
+.idea/**/usage.statistics.xml
+.idea/**/dictionaries
+.idea/**/shelf
+
+# AWS User-specific
+.idea/**/aws.xml
+
+# Generated files
+.idea/**/contentModel.xml
+
+# Sensitive or high-churn files
+.idea/**/dataSources/
+.idea/**/dataSources.ids
+.idea/**/dataSources.local.xml
+.idea/**/sqlDataSources.xml
+.idea/**/dynamic.xml
+.idea/**/uiDesigner.xml
+.idea/**/dbnavigator.xml
+
+# Gradle
+.idea/**/gradle.xml
+.idea/**/libraries
+
+# Gradle and Maven with auto-import
+# When using Gradle or Maven with auto-import, you should exclude module files,
+# since they will be recreated, and may cause churn.  Uncomment if using
+# auto-import.
+# .idea/artifacts
+# .idea/compiler.xml
+# .idea/jarRepositories.xml
+# .idea/modules.xml
+# .idea/*.iml
+# .idea/modules
+# *.iml
+# *.ipr
+
+# CMake
+cmake-build-*/
+
+# Mongo Explorer plugin
+.idea/**/mongoSettings.xml
+
+# File-based project format
+*.iws
+
+# IntelliJ
+out/
+
+# mpeltonen/sbt-idea plugin
+.idea_modules/
+
+# JIRA plugin
+atlassian-ide-plugin.xml
+
+# Cursive Clojure plugin
+.idea/replstate.xml
+
+# SonarLint plugin
+.idea/sonarlint/
+
+# Crashlytics plugin (for Android Studio and IntelliJ)
+com_crashlytics_export_strings.xml
+crashlytics.properties
+crashlytics-build.properties
+fabric.properties
+
+# Editor-based Rest Client
+.idea/httpRequests
+
+# Android studio 3.1+ serialized cache file
+.idea/caches/build_file_checksums.ser
+
+### CMake template
+CMakeLists.txt.user
+CMakeCache.txt
+CMakeFiles
+CMakeScripts
+Testing
+Makefile
+cmake_install.cmake
+install_manifest.txt
+compile_commands.json
+CTestTestfile.cmake
+_deps
+
+### Python template
+# Byte-compiled / optimized / DLL files
+__pycache__/
+*.py[cod]
+*$py.class
+
+# C extensions
+*.so
+
+# Distribution / packaging
+.Python
+build/
+develop-eggs/
+dist/
+downloads/
+eggs/
+.eggs/
+lib/
+lib64/
+parts/
+sdist/
+var/
+wheels/
+share/python-wheels/
+*.egg-info/
+.installed.cfg
+*.egg
+MANIFEST
+
+# PyInstaller
+#  Usually these files are written by a python script from a template
+#  before PyInstaller builds the exe, so as to inject date/other infos into it.
+*.manifest
+*.spec
+
+# Installer logs
+pip-log.txt
+pip-delete-this-directory.txt
+
+# Unit test / coverage reports
+htmlcov/
+.tox/
+.nox/
+.coverage
+.coverage.*
+.cache
+nosetests.xml
+coverage.xml
+*.cover
+*.py,cover
+.hypothesis/
+.pytest_cache/
+cover/
+
+# Translations
+*.mo
+*.pot
+
+# Django stuff:
+*.log
+local_settings.py
+db.sqlite3
+db.sqlite3-journal
+
+# Flask stuff:
+instance/
+.webassets-cache
+
+# Scrapy stuff:
+.scrapy
+
+# Sphinx documentation
+docs/_build/
+
+# PyBuilder
+.pybuilder/
+target/
+
+# Jupyter Notebook
+.ipynb_checkpoints
+
+# IPython
+profile_default/
+ipython_config.py
+
+# pyenv
+#   For a library or package, you might want to ignore these files since the code is
+#   intended to run in multiple environments; otherwise, check them in:
+# .python-version
+
+# pipenv
+#   According to pypa/pipenv#598, it is recommended to include Pipfile.lock in version control.
+#   However, in case of collaboration, if having platform-specific dependencies or dependencies
+#   having no cross-platform support, pipenv may install dependencies that don't work, or not
+#   install all needed dependencies.
+#Pipfile.lock
+
+# poetry
+#   Similar to Pipfile.lock, it is generally recommended to include poetry.lock in version control.
+#   This is especially recommended for binary packages to ensure reproducibility, and is more
+#   commonly ignored for libraries.
+#   https://python-poetry.org/docs/basic-usage/#commit-your-poetrylock-file-to-version-control
+#poetry.lock
+
+# pdm
+#   Similar to Pipfile.lock, it is generally recommended to include pdm.lock in version control.
+#pdm.lock
+#   pdm stores project-wide configurations in .pdm.toml, but it is recommended to not include it
+#   in version control.
+#   https://pdm.fming.dev/#use-with-ide
+.pdm.toml
+
+# PEP 582; used by e.g. github.com/David-OConnor/pyflow and github.com/pdm-project/pdm
+__pypackages__/
+
+# Celery stuff
+celerybeat-schedule
+celerybeat.pid
+
+# SageMath parsed files
+*.sage.py
+
+# Environments
+.env
+.venv
+env/
+venv/
+ENV/
+env.bak/
+venv.bak/
+
+# Spyder project settings
+.spyderproject
+.spyproject
+
+# Rope project settings
+.ropeproject
+
+# mkdocs documentation
+/site
+
+# mypy
+.mypy_cache/
+.dmypy.json
+dmypy.json
+
+# Pyre type checker
+.pyre/
+
+# pytype static type analyzer
+.pytype/
+
+# Cython debug symbols
+cython_debug/
+
+# PyCharm
+#  JetBrains specific template is maintained in a separate JetBrains.gitignore that can
+#  be found at https://github.com/github/gitignore/blob/main/Global/JetBrains.gitignore
+#  and can be added to the global gitignore or merged into this file.  For a more nuclear
+#  option (not recommended) you can uncomment the following to ignore the entire idea folder.
+.idea/
 
-*.[oa]
-cmake-build*/

.gitmodules

@@ -1,12 +0,0 @@
-[submodule "vendor/glfw"]
-	path = vendor/glfw
-	url = https://github.com/glfw/glfw.git
-[submodule "vendor/glm"]
-	path = vendor/glm
-	url = https://github.com/g-truc/glm.git
-[submodule "vendor/toddcox"]
-	path = vendor/toddcox
-	url = https://github.com/JCRaymond/toddcox-faster.git
-[submodule "vendor/yaml-cpp"]
-	path = vendor/yaml-cpp
-	url = git://github.com/jbeder/yaml-cpp.git

.idea/misc.xml

@@ -1,19 +0,0 @@
-<?xml version="1.0" encoding="UTF-8"?>
-<project version="4">
-  <component name="CMakeWorkspace" PROJECT_DIR="$PROJECT_DIR$" />
-  <component name="CidrRootsConfiguration">
-    <sourceRoots>
-      <file path="$PROJECT_DIR$/examples" />
-      <file path="$PROJECT_DIR$/vis/include" />
-      <file path="$PROJECT_DIR$/vis/presets" />
-      <file path="$PROJECT_DIR$/vis/shaders" />
-      <file path="$PROJECT_DIR$/vis/src" />
-    </sourceRoots>
-    <libraryRoots>
-      <file path="$PROJECT_DIR$/vendor" />
-    </libraryRoots>
-  </component>
-  <component name="JavaScriptSettings">
-    <option name="languageLevel" value="ES6" />
-  </component>
-</project>

.idea/toddcox-visualize.iml

@@ -1,2 +0,0 @@
-<?xml version="1.0" encoding="UTF-8"?>
-<module classpath="CMake" type="CPP_MODULE" version="4" />

.idea/vcs.xml

@@ -1,6 +0,0 @@
-<?xml version="1.0" encoding="UTF-8"?>
-<project version="4">
-  <component name="VcsDirectoryMappings">
-    <mapping directory="$PROJECT_DIR$" vcs="Git" />
-  </component>
-</project>

CMakeLists.txt

@@ -1,15 +1,14 @@
 cmake_minimum_required(VERSION 3.10)
 project(toddcox-faster)
 
-set(CMAKE_CXX_STANDARD 17)
+set(CMAKE_CXX_STANDARD 20)
 
-add_subdirectory(vendor/toddcox)
+file(GLOB externals ext/*.cmake)
+foreach (external ${externals})
+    include(${external})
+endforeach ()
 
-add_subdirectory(vendor/yaml-cpp)
-add_subdirectory(vendor/glad)
-add_subdirectory(vendor/glfw)
-add_subdirectory(vendor/glm)
+include(cmake/embed.cmake)
 
+add_subdirectory(tc)
 add_subdirectory(vis)
-
-add_subdirectory(examples)

cmake/embed.cmake

@@ -0,0 +1,86 @@
+find_program(EMBED_LD ${CMAKE_LINKER})
+find_program(EMBED_OBJCOPY ${CMAKE_OBJCOPY})
+
+function(_generate_embed_source EMBED_NAME)
+  set(options)
+  set(oneValueArgs SOURCE HEADER)
+  set(multiValueArgs FILES)
+  cmake_parse_arguments(PARSE "${options}" "${oneValueArgs}" "${multiValueArgs}" ${ARGN})
+
+  set(VIEW_DECLARATIONS)
+  set(VOID_DECLARATIONS)
+  set(DEFINITIONS)
+
+  foreach (FILE IN LISTS PARSE_FILES)
+    get_filename_component(FILE_NAME "${FILE}" NAME)
+
+    string(MAKE_C_IDENTIFIER "_binary_${FILE}" SYMBOL)
+    string(MAKE_C_IDENTIFIER "${FILE_NAME}" IDENTIFIER)
+
+    string(APPEND VIEW_DECLARATIONS
+      "  /// ${FILE} Text Contents\n"
+      "  extern std::string_view const ${IDENTIFIER};\n")
+
+    string(APPEND VOID_DECLARATIONS
+      "  /// ${FILE} Binary Contents\n"
+      "  extern void* const ${IDENTIFIER};\n")
+
+    string(APPEND DEFINITIONS
+      "// ${IDENTIFIER} (${FILE})\n"
+      "extern \"C\" const char ${SYMBOL}_start[], ${SYMBOL}_end[];\n"
+      "std::string_view const ${EMBED_NAME}::${IDENTIFIER}(${SYMBOL}_start, ${SYMBOL}_end);\n"
+      "void* const ${EMBED_NAME}::bin::${IDENTIFIER} = (void *) ${SYMBOL}_start;\n\n")
+  endforeach ()
+
+  file(WRITE "${PARSE_HEADER}"
+    "#pragma once\n"
+    "#include <string>\n\n"
+    "namespace ${EMBED_NAME} {\n${VIEW_DECLARATIONS}}\n\n"
+    "namespace ${EMBED_NAME}::bin {\n${VOID_DECLARATIONS}}\n")
+
+  file(WRITE "${PARSE_SOURCE}"
+    "#include <${EMBED_NAME}.hpp>\n\n"
+    "${DEFINITIONS}")
+endfunction()
+
+function(_embed_file OUTPUT_OBJECT FILE)
+  set(OBJECT "${CMAKE_CURRENT_BINARY_DIR}/${FILE}.o")
+
+  set(${OUTPUT_OBJECT} ${OBJECT} PARENT_SCOPE)
+
+  add_custom_command(
+    COMMENT "Embedding ${FILE} in ${OBJECT}"
+    OUTPUT "${FILE}.o" DEPENDS "${FILE}"
+    WORKING_DIRECTORY "${CMAKE_CURRENT_SOURCE_DIR}"
+    COMMAND ${CMAKE_COMMAND} -E make_directory ${CMAKE_CURRENT_BINARY_DIR}
+    COMMAND ${EMBED_LD} -r -o "${OBJECT}" --format=binary "${FILE}"
+    COMMAND ${EMBED_OBJCOPY} --rename-section .data=.rodata,alloc,load,readonly,data,contents "${OBJECT}"
+    VERBATIM
+  )
+endfunction()
+
+function(add_embed_library EMBED_NAME)
+  set(FILES ${ARGN})
+
+  set(EMBED_ROOT ${CMAKE_CURRENT_BINARY_DIR}/_embed/${EMBED_NAME})
+  set(EMBED_SOURCE "${EMBED_ROOT}/${EMBED_NAME}.cpp")
+  set(EMBED_INCLUDE "${EMBED_ROOT}/include")
+  set(EMBED_HEADER "${EMBED_INCLUDE}/${EMBED_NAME}.hpp")
+
+  set(OBJECTS)
+  foreach (FILE ${ARGN})
+    _embed_file(OBJECT ${FILE})
+    list(APPEND OBJECTS ${OBJECT})
+  endforeach ()
+
+  message(STATUS "Generating embedding library ${EMBED_NAME}")
+  _generate_embed_source(
+    ${EMBED_NAME}
+    SOURCE ${EMBED_SOURCE}
+    HEADER ${EMBED_HEADER}
+    FILES ${FILES})
+
+  add_library(${EMBED_NAME} STATIC ${OBJECTS} "${EMBED_SOURCE}")
+  target_include_directories(${EMBED_NAME} PUBLIC "${EMBED_INCLUDE}")
+  set_target_properties(${EMBED_NAME} PROPERTIES POSITION_INDEPENDENT_CODE ON)
+endfunction()

entt-question.md

@@ -0,0 +1,71 @@
+I have some "composite" objects that I'm not sure how best to map into an ECS. Each composite object has multiple "parts" with various properties. But the parts of a given composite object are intrinsically linked: they transform together and use the same buffers during rendering.
+
+Note that I'm still pretty early in adopting `entt` for this project, so I'm not particularly tied to any of the architecture here. The only constraint is that, for optimization reasons, I need each composite object to get _one_ buffer on GPU, and each part owns some interval within its parent's buffer. This nested scheme seems the most straightforward way to do that.
+
+---
+
+My first attempt has only one entity for each composite object, with a `Parts` component that holds an `std::vector<Part>`. 
+
+Rendering code looks something like this. In reality I save the draw commands in a separate component and re-use those on each frame.
+
+```c++
+auto view = registry.view<Parts>();
+
+for (auto [entity, parts]: view.each()) {
+  std::vector<Command> commands;
+
+  for (auto part: parts.parts) {
+    // assemble draw command
+    commands.emplace_back(...);
+  }
+
+  // bind buffers for entity
+  // issue draw command
+}
+```
+
+However it seems appropriate to give each part its own entity. My attempt for this was to switch to `std::vector<entt::entity>` where each entity has a `Part` as a component.
+
+Rendering code then looks like:
+
+```c++
+auto view = registry.view<Parts>();
+
+for (auto [entity, parts]: view.each()) {
+  std::vector<Command> commands;
+
+  for (auto part_entity: parts.parts) {
+    auto &part = registry.get<Part>(part_entity);
+    // assemble draw command
+    commands.emplace_back(...);
+  }
+
+  // bind buffers for entity
+  // issue draw command with shared buffer
+}
+```
+
+---
+
+In my last attempt, I invert this so each `Part` holds only parent entity. The `Parts` component is removed. I can then use a `view<Part>` directly, but it's trickier to organize the draw calls since they must be grouped by the parent entity itself.
+
+```c++
+auto view = registry.view<Part>();
+
+for (auto [entity, part]: view.each()) {
+  // assemble draw command
+  auto &commands = registry.get<Commands>(part.parent);
+  commands.emplace_back(...);
+}
+
+for (auto [entity, commands]: view
+for (auto [entity, commands]: commands) {
+   // bind buffers for entity
+   auto &commands = registry.get<Commands>(entity);
+   // issue draw commands
+}
+```
+
+---
+
+Is there any more sensible approach to this kind of many-to-one entity composition? Or some way to declare _multiple_ `Part` components on the same entity? It seems like I almost need some group-by kind of operation, but I don't see such on the registry/view docs.

examples/CMakeLists.txt

@@ -1 +0,0 @@
-# there's no examples right now.

ext/eigen.cmake

@@ -0,0 +1,14 @@
+include(FetchContent)
+
+FetchContent_Declare(
+    eigen
+    GIT_REPOSITORY https://gitlab.com/libeigen/eigen.git
+    GIT_TAG 3.4
+    GIT_PROGRESS TRUE
+)
+
+set(EIGEN_BUILD_DOC OFF CACHE INTERNAL "")
+set(BUILD_TESTING OFF CACHE INTERNAL "")
+set(EIGEN_BUILD_PKGCONFIG OFF CACHE INTERNAL "")
+
+FetchContent_MakeAvailable(eigen)

ext/entt.cmake

@@ -0,0 +1,10 @@
+include(FetchContent)
+
+FetchContent_Declare(
+    entt
+    GIT_REPOSITORY https://github.com/skypjack/entt.git
+    GIT_TAG v3.11.1
+    GIT_PROGRESS TRUE
+)
+
+FetchContent_MakeAvailable(entt)

ext/fmt.cmake

@@ -0,0 +1,10 @@
+include(FetchContent)
+
+FetchContent_Declare(
+    fmt
+    GIT_REPOSITORY https://github.com/fmtlib/fmt.git
+    GIT_TAG 9.1.0
+    GIT_PROGRESS TRUE
+)
+
+FetchContent_MakeAvailable(fmt)

ext/glad.cmake

@@ -0,0 +1,14 @@
+include(FetchContent)
+
+FetchContent_Declare(
+    glad
+    GIT_REPOSITORY https://github.com/Dav1dde/glad.git
+    GIT_TAG v0.1.36
+    GIT_PROGRESS TRUE
+)
+
+set(GLAD_PROFILE "core" CACHE INTERNAL "OpenGL profile")
+set(GLAD_API "gl=4.6" CACHE INTERNAL "API type/version pairs, like \"gl=3.2,gles=\", no version means latest")
+set(GLAD_GENERATOR "c" CACHE INTERNAL "Language to generate the binding for")
+
+FetchContent_MakeAvailable(glad)

ext/glfw.cmake

@@ -0,0 +1,15 @@
+include(FetchContent)
+
+FetchContent_Declare(
+    glfw
+    GIT_REPOSITORY https://github.com/glfw/glfw.git
+    GIT_TAG 3.3.8
+    GIT_PROGRESS TRUE
+)
+
+set(GLFW_BUILD_EXAMPLES OFF CACHE INTERNAL "")
+set(GLFW_BUILD_TESTS OFF CACHE INTERNAL "")
+set(GLFW_BUILD_DOCS OFF CACHE INTERNAL "")
+set(GLFW_INSTALL OFF CACHE INTERNAL "")
+
+FetchContent_MakeAvailable(glfw)

ext/gtest.cmake

@@ -0,0 +1,12 @@
+include(FetchContent)
+
+FetchContent_Declare(
+    gtest
+    GIT_REPOSITORY https://github.com/google/googletest.git
+    GIT_TAG release-1.12.1
+    GIT_PROGRESS TRUE
+)
+
+set(gtest_force_shared_crt ON CACHE BOOL "" FORCE)
+
+FetchContent_MakeAvailable(gtest)

ext/imgui.cmake

@@ -0,0 +1,31 @@
+include(FetchContent)
+
+FetchContent_Declare(
+    imgui
+    GIT_REPOSITORY https://github.com/ocornut/imgui
+    GIT_TAG v1.86
+    GIT_PROGRESS TRUE
+)
+
+FetchContent_GetProperties(imgui)
+if (NOT ${imgui}_POPULATED)
+    FetchContent_MakeAvailable(imgui)
+
+    add_library(
+        imgui
+        ${imgui_SOURCE_DIR}/backends/imgui_impl_glfw.cpp
+        ${imgui_SOURCE_DIR}/backends/imgui_impl_glfw.h
+        ${imgui_SOURCE_DIR}/backends/imgui_impl_opengl3.cpp
+        ${imgui_SOURCE_DIR}/backends/imgui_impl_opengl3.h
+        ${imgui_SOURCE_DIR}/imgui.cpp
+        ${imgui_SOURCE_DIR}/imgui.h
+        ${imgui_SOURCE_DIR}/imgui_demo.cpp
+        ${imgui_SOURCE_DIR}/imgui_draw.cpp
+        ${imgui_SOURCE_DIR}/imgui_tables.cpp
+        ${imgui_SOURCE_DIR}/imgui_widgets.cpp
+    )
+
+    target_include_directories(imgui PUBLIC ${imgui_SOURCE_DIR})
+    target_link_libraries(imgui PRIVATE glfw)
+endif ()
+

ext/peglib.cmake

@@ -0,0 +1,20 @@
+include(FetchContent)
+
+FetchContent_Declare(
+    peglib
+    GIT_REPOSITORY https://github.com/yhirose/cpp-peglib
+    GIT_TAG v1.8.2
+    GIT_PROGRESS TRUE
+)
+
+set(PEGLIB_BUILD_TESTS OFF CACHE INTERNAL "")
+
+FetchContent_GetProperties(peglib)
+if(NOT ${peglib}_POPULATED)
+    FetchContent_Populate(peglib)
+    find_package(Threads)
+    add_library(peglib INTERFACE ${peglib_SOURCE_DIR}/peglib.h)
+    target_include_directories(peglib INTERFACE ${peglib_SOURCE_DIR})
+    target_link_libraries(peglib INTERFACE Threads::Threads)
+    add_library(peglib::peglib ALIAS peglib)
+endif()

ext/yaml-cpp.cmake

@@ -0,0 +1,15 @@
+include(FetchContent)
+
+FetchContent_Declare(
+    yaml-cpp
+    GIT_REPOSITORY https://github.com/jbeder/yaml-cpp.git
+    GIT_TAG yaml-cpp-0.7.0
+    GIT_PROGRESS TRUE
+)
+
+set(GLFW_BUILD_EXAMPLES OFF CACHE INTERNAL "")
+set(GLFW_BUILD_TESTS OFF CACHE INTERNAL "")
+set(GLFW_BUILD_DOCS OFF CACHE INTERNAL "")
+set(GLFW_INSTALL OFF CACHE INTERNAL "")
+
+FetchContent_MakeAvailable(yaml-cpp)

tc/CMakeLists.txt

@@ -0,0 +1,17 @@
+add_library(tc
+    include/tc/core.hpp
+    include/tc/groups.hpp
+
+    src/cosets.cpp
+    src/group.cpp
+    src/groups.cpp
+    src/lang.cpp
+    src/solve.cpp
+    )
+target_link_libraries(tc peglib fmt)
+target_include_directories(tc PUBLIC include)
+
+add_library(tc::tc ALIAS tc)
+
+add_subdirectory(test)
+add_subdirectory(bench)

tc/bench/CMakeLists.txt

@@ -0,0 +1,5 @@
+add_executable(benchmark benchmark.cpp)
+target_link_libraries(benchmark PUBLIC tc fmt::fmt)
+
+add_executable(named named.cpp)
+target_link_libraries(named PUBLIC tc fmt::fmt)

tc/bench/benchmark.cpp

@@ -0,0 +1,138 @@
+#include <ctime>
+#include <string>
+#include <vector>
+
+#include <fmt/core.h>
+#include <fmt/ranges.h>
+
+#include <tc/core.hpp>
+#include <tc/groups.hpp>
+
+void bench(
+    const std::string &group_expr,
+    const std::string &symbol,
+    const std::vector<size_t> &gens,
+    const size_t bound = SIZE_MAX
+) {
+    tc::Group<> group = tc::coxeter(symbol);
+    
+    std::clock_t s = std::clock();
+    tc::Cosets<> cosets = group.solve(gens, bound);
+    std::clock_t e = std::clock();
+
+    auto time = (double) (e - s) / CLOCKS_PER_SEC;
+    size_t order = cosets.order();
+    auto cos_s = (size_t) (order / time);
+
+    bool complete = cosets.complete();
+
+    std::string name = fmt::format("{}/{}", group_expr, gens);
+    std::string row = fmt::format(
+        "{:>24},{:>10},{:>6},{:>8.3f}s,{:>10L}",
+        name, order, complete, time, cos_s
+    );
+    fmt::print("{}\n", row);
+}
+
+int main(int argc, char *argv[]) {
+    std::vector<std::string> args(argv + 1, argv + argc);
+
+    fmt::print("{:>24},{:>10},{:>6},{:>9},{:>10}\n", "NAME", "ORDER", "COMPL", "TIME", "COS/S");
+
+    // Finite Groups
+    
+    // A_n: 3 * `n-1`           ; n >= 1
+    bench("A_5", "3 * 4", {});
+    bench("A_6", "3 * 5", {});
+    bench("A_7", "3 * 6", {});
+    bench("A_8", "3 * 7", {});
+    // B_n: 4 3 * `n-2`         ; n >= 2
+    bench("B_5", "4 3 * 3", {});
+    bench("B_6", "4 3 * 4", {});
+    bench("B_7", "4 3 * 5", {});
+    bench("B_8", "4 3 * 6", {});
+    // D_n: 3 * [1 1 `n-3`]     ; n >= 4
+    bench("D_5", "3 * [1 1 2]", {});
+    bench("D_6", "3 * [1 1 3]", {});
+    bench("D_7", "3 * [1 1 4]", {});
+    bench("D_8", "3 * [1 1 5]", {});
+    // E_n: 3 * [1 2 `n-4`]     ; n >= 6
+    bench("E_6", "3 * [1 2 2]", {});
+    bench("E_7", "3 * [1 2 3]", {});
+//    bench("E_8", "3 * [1 2 4]", {}); // too big
+    // H_n: 5 3 * `n-2`         ; n >= 2
+    bench("H_3", "5 3 * 0", {});
+    bench("H_4", "5 3 * 1", {});
+    bench("H_5", "5 3 * 2", {});
+    // grid: `p` `q`            ; 2(p+q) > pq
+    // triangle: `p` `q` `r`    ; 1/p + 1/q + 1/r > 1
+
+    // Special Finite Groups
+    bench("F_4", "3 4 3", {});
+    bench("G_2", "6", {});
+    // I_2(p): `p`              ; p >= 2
+    bench("I_2(100)", "100", {});
+    bench("I_2(1000)", "1000", {});
+    // "Torus": `p` 2 `q`       ; p, q >= 2
+    bench("T(100)", "100 2 100", {});
+    bench("T(1000)", "1000 2 1000", {});
+
+    // Affine Groups
+    
+    // ~A_n: {3 * `n+1`}
+    bench("~A_5", "{3 * 6}", {}, 10'000'000);
+    bench("~A_6", "{3 * 7}", {}, 10'000'000);
+    bench("~A_7", "{3 * 8}", {}, 10'000'000);
+    bench("~A_8", "{3 * 9}", {}, 10'000'000);
+    // ~B_n: 4 3 * `n-3` 3 * [1 1]
+    bench("~B_5", "4 3 * 2 3 * [1 1]", {}, 10'000'000);
+    bench("~B_6", "4 3 * 3 3 * [1 1]", {}, 10'000'000);
+    bench("~B_7", "4 3 * 4 3 * [1 1]", {}, 10'000'000);
+    bench("~B_8", "4 3 * 5 3 * [1 1]", {}, 10'000'000);
+    // ~B_n: 4 3 * `n-2` 4
+    bench("~C_5", "4 3 * 3 4", {}, 10'000'000);
+    bench("~C_6", "4 3 * 4 4", {}, 10'000'000);
+    bench("~C_7", "4 3 * 5 4", {}, 10'000'000);
+    bench("~C_8", "4 3 * 6 4", {}, 10'000'000);
+    // ~D_n: 3 * [1 1] 3 * `n-4` 3 * [1 1]
+    bench("~D_5", "3 * [1 1] 3 * 1 3 * [1 1]", {}, 10'000'000);
+    bench("~D_6", "3 * [1 1] 3 * 2 3 * [1 1]", {}, 10'000'000);
+    bench("~D_7", "3 * [1 1] 3 * 3 3 * [1 1]", {}, 10'000'000);
+    bench("~D_8", "3 * [1 1] 3 * 4 3 * [1 1]", {}, 10'000'000);
+    // grid: `p` `q`            ; 2(p+q) = pq
+    // triangle: `p` `q` `r`    ; 1/p + 1/q + 1/r = 1
+    
+    // Special Affine Groups
+    bench("~I_1", "-", {}, 10'000'000);
+    bench("~E_6", "3 * [2 2 2]", {}, 10'000'000);
+    bench("~E_7", "3 * [1 3 3]", {}, 10'000'000);
+    bench("~E_8", "3 * [1 2 5]", {}, 10'000'000);
+//    bench("E_9",  "3 * [1 2 5]", {}, 10'000'000);  // ~E_8 == E_9
+    bench("~F_4", "3 4 3 3", {}, 10'000'000);
+    bench("~G_2", "6 3", {}, 10'000'000);
+
+    // Hyperbolic Groups
+    // grid: `p` `q`            ; 2(p+q) < pq
+    // triangle: `p` `q` `r`    ; 1/p + 1/q + 1/r < 1
+    
+    // Special Hyperbolic Groups
+    bench("-BH_3", "4 3 5", {}, 10'000'000);
+    bench("-K_3", "5 3 5", {}, 10'000'000);
+    bench("-J_3", "3 5 3", {}, 10'000'000);
+//    bench("~H_3", "3 5 3", {}, 10'000'000);  // -J_3 == ~H_3
+    bench("-DH_3", "5 3 * [1 1]", {}, 10'000'000);
+    bench("^AB_3", "{3 3 3 4}", {}, 10'000'000);
+    bench("^AH_3", "{3 3 3 5}", {}, 10'000'000);
+    bench("^BB_3", "{3 4 3 4}", {}, 10'000'000);
+    bench("^BH_3", "{3 4 3 5}", {}, 10'000'000);
+    bench("^HH_3", "{3 5 3 5}", {}, 10'000'000);
+    bench("-H_4", "5 3 3 3", {}, 10'000'000);
+//    bench("~H_4", "5 3 3 3", {}, 10'000'000);  // -H_4 == ~H_4 == H_5 
+//    bench("H_5", "5 3 3 3", {}, 10'000'000);
+    bench("-BH_4", "4 3 3 5", {}, 10'000'000);
+    bench("-K_4", "5 3 3 5", {}, 10'000'000);
+    bench("-DH_4", "5 3 3 * [1 1]", {}, 10'000'000);
+    bench("^AF_4", "{3 3 3 3 4}", {}, 10'000'000);
+    
+    return EXIT_SUCCESS;
+}

tc/bench/named.cpp

@@ -0,0 +1,149 @@
+#include <string>
+#include <vector>
+
+#include <fmt/core.h>
+#include <fmt/ranges.h>
+
+#include <tc/core.hpp>
+
+using namespace std::string_literals;
+
+template<typename Gen_>
+struct fmt::formatter<tc::Group<Gen_>> {
+    template<typename ParseContext>
+    constexpr auto parse(ParseContext &ctx) {
+        return ctx.begin();
+    }
+
+    template<typename FormatContext>
+    constexpr auto format(tc::Group<Gen_> const &g, FormatContext &ctx) {
+        auto gens = g.gens();
+
+        fmt::format_to(ctx.out(), "  | ");
+        for (const auto &gen: gens) {
+            fmt::format_to(ctx.out(), "{} ", gen);
+        }
+        fmt::format_to(ctx.out(), "\n");
+
+        for (int i = 0; i < gens.size(); ++i) {
+            auto u = gens[i];
+
+            fmt::format_to(ctx.out(), "{} | ", u);
+
+            for (int j = 0; j < gens.size(); ++j) {
+                auto v = gens[j];
+
+                if (i <= j) {
+                    fmt::format_to(ctx.out(), "{} ", g.get(u, v));
+                } else {
+                    fmt::format_to(ctx.out(), "  ");
+                }
+            }
+
+            fmt::format_to(ctx.out(), "\n");
+        }
+
+        return ctx.out();
+    }
+};
+
+template<typename Gen_>
+struct fmt::formatter<tc::Cosets<Gen_>> {
+    template<typename ParseContext>
+    constexpr auto parse(ParseContext &ctx) {
+        return ctx.begin();
+    }
+
+    template<typename FormatContext>
+    constexpr auto format(tc::Cosets<Gen_> const &c, FormatContext &ctx) {
+        auto gens = c.gens();
+
+        auto width = (size_t) (std::log10(c.order() + 1));
+
+        fmt::format_to(ctx.out(), "{:>{}} | ", "-", width);
+        for (const auto &gen: gens) {
+            fmt::format_to(ctx.out(), "{:>{}} ", gen, width);
+        }
+        fmt::format_to(ctx.out(), "\n");
+
+        for (size_t cos = 0; cos < c.order(); ++cos) {
+            fmt::format_to(ctx.out(), "{:>{}} | ", cos, width);
+
+            for (const auto &gen: gens) {
+                auto target = c.get(cos, gen);
+
+                fmt::format_to(ctx.out(), "{:>{}} ", target, width);
+            }
+
+            fmt::format_to(ctx.out(), "\n");
+        }
+
+        return ctx.out();
+    }
+};
+
+int main() {
+    tc::Group<char> group(4, {'r', 'g', 'b', 'y'});
+
+    group.set('r', 'g', 5);
+    group.set('g', 'b', 3);
+    group.set('b', 'y', 3);
+
+    fmt::print("{}\n", group);
+
+//    auto res = group.sub({'r', 'g', 'b', 'y'}).solve();
+    auto res = group.sub({'r', 'g'}).solve();
+    fmt::print("rank: {}, order: {}\n{}\n",
+               res.rank(), res.order(), res);
+
+//    tc::Path<> path(res);
+//    auto gens = group.gens();
+//    std::vector<std::string> words(path.order());
+//    path.walk(
+//        "-"s,
+//        [&](const std::string &cos, size_t gen) -> std::string {
+//            return cos + gens[gen];
+//        },
+//        words.begin());
+
+    tc::Path<char> typed(res);
+    std::vector<std::string> words(typed.order());
+    typed.walk("-"s, std::plus<>(), words.begin());
+
+    fmt::print("words: {}\n", words);
+    fmt::print("size: {}\n", words.size());
+
+//    path.walk(
+//        words.begin(),
+//        "-",
+//        [&gens](const std::string &cos, size_t gen) {
+//            return cos + gens[gen];
+//        }
+//    );
+//    for (const auto &word: words) {
+//        fmt::print("'{}'\n", word);
+//    }
+
+//    auto sub = group.sub({'r', 'g', 'y'});
+//    fmt::print("{}\n", sub);
+//
+//    auto res = sub.solve({});
+//    fmt::print("res order: {}\n", res.order());
+//
+//    auto cos = sub.solve({'r', 'y'});
+
+//    fmt::print("order: {}\n", cos.order());
+//
+//    fmt::print("{}\n", cos);
+
+//    tc::Group<u_int8_t> group(4, {0, 1, 2, 3});
+//
+//    group.set(0, 1, 5);
+//    group.set(1, 2, 4);
+//    group.set(2, 3, 3);
+//
+//    show(group);
+//
+//    auto sub = group.sub({3, 2, 0, 1});
+//    show(sub);
+}

tc/include/tc/core.hpp

@@ -0,0 +1,194 @@
+#pragma once
+
+#include <vector> //todo clean up includes. lots of duplicate cstdint, cassert.
+#include <cstdint>
+#include <cassert>
+#include <tuple>
+
+#include <limits>
+
+namespace tc {
+    using Mult = std::uint16_t;
+    constexpr Mult FREE = 0;
+
+    /**
+     * @brief Complete representation of a quotient group. Describes the action of each generator on each coset.
+     */
+    struct Cosets;
+
+    /** 
+     * @brief Manage the presentation of a Coxeter group and enforce constraints
+     * on the multiplicities of its relations.
+     * <ul>
+     *   <li>
+     *     <code>m_ij = 1</code> iff <code>i != j</code>
+     *   </li>
+     *   <li>
+     *     <code>m_ij = m_ji</code>
+     *   </li>
+     *   <li>
+     *     If <code>m_ij == inf</code> (<code>tc::FREE</code>) then no relation is imposed.
+     *   </li>
+     * </ul>
+     * @see
+     * <a href="https://en.wikipedia.org/wiki/Coxeter_group#Definition">Coxeter Group (Wikipedia)</a>
+     */
+    struct Group;
+
+    /**
+     * @brief Support generating values given a Cosets and transformation callback.
+     * @tparam Gen_ 
+     */
+    template<typename Gen_=void>
+    struct Path;
+
+    struct Cosets {
+        static constexpr size_t UNSET = std::numeric_limits<size_t>::max();
+
+    private:
+        size_t _rank;
+        size_t _order;
+        bool _complete;
+        std::vector<size_t> _data;
+
+    public:
+        Cosets(Cosets const &) = default;
+
+        Cosets(Cosets &&) noexcept = default;
+
+        ~Cosets() = default;
+
+        void set(size_t coset, size_t gen, size_t target);
+
+        [[nodiscard]] size_t get(size_t coset, size_t gen) const;
+
+        [[nodiscard]] bool isset(size_t coset, size_t gen) const;
+
+        [[nodiscard]] size_t rank() const;
+
+        [[nodiscard]] size_t order() const;
+
+        [[nodiscard]] bool complete() const;
+
+        [[nodiscard]] size_t size() const;
+
+        friend Group;  // only constructible via Group<>::solve
+
+    private:
+        explicit Cosets(size_t rank);
+
+        void add_row();
+
+        void set(size_t idx, size_t target);
+
+        [[nodiscard]] size_t get(size_t idx) const;
+
+        [[nodiscard]] bool isset(size_t idx) const;
+    };
+
+    struct Group {
+        using Rel = std::tuple<size_t, size_t, Mult>;
+
+    private:
+        size_t _rank;
+        std::vector<size_t> _mults;
+        std::vector<size_t> _gens;
+
+    public:
+        Group(Group const &) = default;
+
+        Group(Group &&) noexcept = default;
+
+        ~Group() = default;
+
+        explicit Group(size_t rank);
+
+        void set(size_t, size_t, Mult);
+
+        [[nodiscard]] Mult get(size_t, size_t) const;
+
+        [[nodiscard]] size_t rank() const;
+
+        [[nodiscard]] std::vector<size_t> gens() const;
+
+        [[nodiscard]] Group sub(std::vector<size_t> const &gens) const;
+
+        [[nodiscard]] std::vector<Group> subs(size_t rank) const;
+
+        [[nodiscard]] Cosets solve(std::vector<size_t> const &idxs = {}, size_t bound = SIZE_MAX) const;
+    };
+
+    template<>
+    struct Path<> {
+        using Action = std::tuple<size_t, size_t>;  // coset, gen
+
+    private:
+        std::vector<Action> _data;
+
+    public:
+        // todo be smarter about move semantics
+
+        explicit Path(Cosets const &cosets) : _data() {
+            _data.resize(cosets.order());
+
+            std::vector<bool> complete(cosets.order(), false);
+
+            complete[0] = true;
+
+            for (size_t cos = 0; cos < cosets.order(); ++cos) {
+                for (size_t gen = 0; gen < cosets.rank(); ++gen) {
+                    size_t tgt = cosets.get(cos, gen);
+                    if (!complete[tgt]) {
+                        _data[tgt] = {cos, gen};
+                        complete[tgt] = true;
+                    }
+                }
+            }
+        }
+
+        template<typename Elem, typename BinaryOp, std::random_access_iterator It>
+        void walk(Elem const &start, BinaryOp op, It out) {
+            // todo how to work for std::back_insert_iterator?
+            // todo how to work for non-default-constructible values?
+
+            out[0] = start;
+
+            for (int tgt = 1; tgt < _data.size(); ++tgt) {
+                auto [cos, gen] = _data[tgt];
+                out[tgt] = op(out[cos], gen);
+            }
+        }
+
+        template<typename Elem, typename BinaryOp, std::random_access_iterator It, std::random_access_iterator Gens>
+        void walk(Elem const &start, BinaryOp op, It out, Gens const &gens) {
+            out[0] = start;
+
+            for (int tgt = 1; tgt < _data.size(); ++tgt) {
+                auto [cos, gen] = _data[tgt];
+                out[tgt] = op(out[cos], gens[gen]);
+            }
+        }
+
+        [[nodiscard]] size_t order() const {
+            return _data.size();
+        }
+    };
+
+    template<typename Gen_>
+    struct Path : public Path<> {
+        using Gen = Gen_;
+
+    private:
+        std::vector<Gen> _index;
+
+    public:
+        template<typename T>
+        Path(Cosets const &cosets, T &&gens)
+            : Path<>(cosets), _index(std::forward<T>(gens).begin(), std::forward<T>(gens).end()) {}
+
+        template<typename Elem, typename BinaryOp, std::random_access_iterator It>
+        void walk(Elem const &start, BinaryOp op, It out) {
+            Path<>::walk(start, op, out, _index.begin());
+        }
+    };
+}

tc/include/tc/groups.hpp

@@ -0,0 +1,24 @@
+#pragma once
+
+#include <vector>
+#include <string>
+
+#include <tc/core.hpp>
+
+namespace tc {
+    /**
+     * Construct a group from a (simplified) Schlafli Symbol of the form [a, b, ..., c]
+     * @param mults: The sequence of multiplicites between adjacent generators.
+     */
+    Group schlafli(const std::vector<unsigned int> &mults);
+
+    Group coxeter(const std::string &symbol);
+
+    Group vcoxeter(const std::string &symbol, const std::vector<unsigned int> &values);
+
+    template<typename ...Args>
+    Group coxeter(const std::string &symbol, const Args &... args) {
+        std::vector<unsigned int> values = {{args...}};
+        return vcoxeter(symbol, values);
+    }
+}

tc/src/cosets.cpp

@@ -0,0 +1,55 @@
+#include <tc/core.hpp>
+
+namespace tc {
+    Cosets::Cosets(size_t rank)
+        : _rank(rank), _order(0), _complete(false), _data() {}
+
+    void Cosets::set(size_t coset, size_t gen, size_t target) {
+        set(coset * rank() + gen, target);
+    }
+
+    [[nodiscard]] size_t Cosets::get(size_t coset, size_t gen) const {
+        return get(coset * rank() + gen);
+    }
+
+    [[nodiscard]] bool Cosets::isset(size_t coset, size_t gen) const {
+        return isset(coset * rank() + gen);
+    }
+
+    [[nodiscard]] size_t Cosets::rank() const {
+        return _rank;
+    }
+
+    [[nodiscard]] size_t Cosets::order() const {
+        return _order;
+    }
+
+    [[nodiscard]] bool Cosets::complete() const {
+        return _complete;
+    }
+
+    [[nodiscard]] size_t Cosets::size() const {
+        return _data.size();
+    }
+
+    void Cosets::add_row() {
+        _data.resize(_data.size() + rank(), UNSET);
+        _order++;
+    }
+
+    void Cosets::set(size_t idx, size_t target) {
+        size_t coset = idx / rank();
+        size_t gen = idx % rank();
+        _data[idx] = target;
+        _data[target * rank() + gen] = coset;
+    }
+
+    [[nodiscard]] size_t Cosets::get(size_t idx) const {
+        return _data[idx];
+    }
+
+    [[nodiscard]] bool Cosets::isset(size_t idx) const {
+        return get(idx) != UNSET;
+    }
+
+}

tc/src/group.cpp

@@ -0,0 +1,72 @@
+#include <tc/core.hpp>
+
+#include <cassert>
+#include <numeric>
+
+namespace tc {
+    Group::Group(size_t rank) : _rank(rank), _mults(_rank * _rank, 2), _gens(_rank) {
+        std::iota(_gens.begin(), _gens.end(), 0);
+
+        for (int idx = 0; idx < rank; ++idx) {
+            set(idx, idx, 1);
+        }
+    }
+
+    void Group::set(size_t u, size_t v, Mult m) {
+        assert(u < rank());
+        assert(v < rank());
+
+        _mults[u * rank() + v] = m;
+        _mults[v * rank() + u] = m;
+    }
+
+    [[nodiscard]] Mult Group::get(size_t u, size_t v) const {
+        assert(u < rank());
+        assert(v < rank());
+
+        return _mults[u * rank() + v];
+    }
+
+    [[nodiscard]] size_t Group::rank() const {
+        return _rank;
+    }
+
+    [[nodiscard]] std::vector<size_t> Group::gens() const {
+        return _gens;
+    }
+
+    [[nodiscard]] Group Group::sub(std::vector<size_t> const &gens) const {
+        Group res(gens.size());
+        res._gens = gens;
+
+        for (int i = 0; i < gens.size(); ++i) {
+            for (int j = i; j < gens.size(); ++j) {
+                res.set(i, j, get(gens[i], gens[j]));
+            }
+        }
+
+        return res;
+    }
+
+    [[nodiscard]] std::vector<Group> Group::subs(size_t rank) const {
+        std::vector<bool> mask(_rank, false);
+        std::fill(mask.begin(), mask.begin() + rank, true);
+
+        std::vector<size_t> sub_gens(rank);
+
+        std::vector<Group> res;
+
+        do {
+            for (int i = 0, j = 0; i < _rank; ++i) {
+                if (mask[i]) sub_gens[j++] = i;
+            }
+            res.push_back(sub(sub_gens));
+        } while (std::next_permutation(
+            mask.begin(),
+            mask.end(),
+            std::greater<>()
+        ));
+
+        return res;
+    }
+}

tc/src/groups.cpp

@@ -0,0 +1,25 @@
+#include <tc/groups.hpp>
+
+#include <fmt/args.h>
+#include <fmt/core.h>
+#include <numeric>
+
+namespace tc {
+    Group schlafli(const std::vector<unsigned int> &mults) {
+        Group res(mults.size() + 1);
+        for (size_t i = 0; i < mults.size(); ++i) {
+            res.set(i, i + 1, mults[i]);
+        }
+        return res;
+    }
+
+    Group vcoxeter(const std::string &symbol, const std::vector<unsigned int> &values) {
+        fmt::dynamic_format_arg_store<fmt::format_context> ds;
+
+        for (const auto &value: values) {
+            ds.push_back(value);
+        }
+
+        return coxeter(fmt::vformat(symbol, ds));
+    }
+}

tc/src/lang.cpp

@@ -0,0 +1,299 @@
+#include <vector>
+#include <stack>
+#include <string>
+#include <cassert>
+
+#include <tc/core.hpp>
+#include <tc/groups.hpp>
+
+#include <peglib.h>
+
+#include <fmt/core.h>
+#include <fmt/ranges.h>
+#include <numeric>
+
+struct Graph {
+    size_t rank{};
+    std::vector<tc::Group::Rel> edges{};
+};
+
+struct Op {
+    enum Code {
+        LINK,
+        PUSH,
+        POP,
+        LOOP,
+        FREE,
+    };
+
+    Code code: 4;
+    unsigned int value: 12;
+
+    explicit Op(Code code, unsigned int value = 0)
+        : code(code), value(value) {}
+};
+
+template<>
+struct fmt::formatter<Op> {
+    template<typename ParseContext>
+    constexpr auto parse(ParseContext &ctx) {
+        return ctx.begin();
+    }
+
+    template<typename FormatContext>
+    auto format(const Op &op, FormatContext &ctx) {
+        switch (op.code) {
+            case Op::LINK:
+                return fmt::format_to(ctx.out(), "link({})",
+                                      (unsigned int) op.value);
+            case Op::PUSH:
+                return fmt::format_to(ctx.out(), "push()");
+            case Op::POP:
+                return fmt::format_to(ctx.out(), "pop()");
+            case Op::LOOP:
+                return fmt::format_to(ctx.out(), "loop()");
+            case Op::FREE:
+                return fmt::format_to(ctx.out(), "free()");
+            default:
+                return fmt::format_to(ctx.out(), "[{}]({})",
+                                      (unsigned int) op.code,
+                                      (unsigned int) op.value);
+        }
+    }
+};
+
+struct Factor {
+    unsigned int mode;
+    std::vector<unsigned int> orders;
+};
+
+struct codegen {
+    std::vector<Op> ops;
+
+    void link(unsigned int order) {
+        ops.emplace_back(Op::LINK, order);
+    }
+
+    void push() {
+        ops.emplace_back(Op::PUSH);
+    }
+
+    void pop() {
+        ops.emplace_back(Op::POP);
+    }
+
+    void loop() {
+        ops.emplace_back(Op::LOOP);
+    }
+
+    void free() {
+        ops.emplace_back(Op::FREE);
+    }
+
+    template<typename It>
+    void insert(It begin, It end) {
+        std::copy(begin, end, std::back_inserter(ops));
+    }
+
+    void insert(const codegen &o) {
+        insert(o.ops.begin(), o.ops.end());
+    }
+};
+
+static const std::string GRAMMAR = R"(
+    root    <- term+
+    term    <- product / op
+    op      <- block / link
+    block   <- '(' root ')' / '{' root '}' / '[' root ']'
+    link    <- int / '-'
+    product <- op '*' factor
+    factor  <- int / '{' int+ '}' / '[' int+ ']'
+    int     <- < [0-9]+ >
+
+    %whitespace <- [ \t\n\r]*
+)";
+
+peg::parser build_parser() {
+    peg::parser parser;
+    parser.set_logger([](size_t line, size_t col, const std::string &msg, const std::string &rule) {
+        fmt::print(stderr, "{}:{} [{}] {}\n", line, col, rule, msg);
+    });
+    auto ok = parser.load_grammar(GRAMMAR);
+    assert(ok);
+
+    parser["int"] = [](const peg::SemanticValues &vs) -> std::any {
+        return vs.token_to_number<unsigned int>();
+    };
+
+    parser["link"] = [](const peg::SemanticValues &vs) -> std::any {
+        codegen cg;
+
+        if (vs.choice() == 0) {
+            auto order = std::any_cast<unsigned int>(vs[0]);
+            cg.link(order);
+        } else {
+            cg.free();
+        }
+
+        return cg;
+    };
+
+    parser["root"] = [](const peg::SemanticValues &vs) -> std::any {
+        codegen cg;
+
+        for (const auto &sub: vs) {
+            cg.insert(std::any_cast<codegen>(sub));
+        }
+
+        return cg;
+    };
+
+    parser["block"] = [](const peg::SemanticValues &vs) -> std::any {
+        if (vs.choice() == 0) return vs[0];
+
+        codegen cg;
+
+        cg.push();
+        cg.insert(std::any_cast<codegen>(vs[0]));
+        if (vs.choice() == 1) cg.loop();
+        cg.pop();
+
+        return cg;
+    };
+
+    parser["factor"] = [](const peg::SemanticValues &vs) -> std::any {
+        return Factor{
+            (unsigned int) vs.choice(),
+            vs.transform<unsigned int>(),
+        };
+    };
+
+    parser["product"] = [](const peg::SemanticValues &vs) -> std::any {
+        auto sub = std::any_cast<codegen>(vs[0]);
+        auto fac = std::any_cast<Factor>(vs[1]);
+
+        codegen cg;
+
+        for (const auto &order: fac.orders) {
+            if (fac.mode == 0) {
+                for (unsigned int i = 0; i < order; ++i) {
+                    cg.insert(sub);
+                }
+            } else {
+                cg.push();
+                for (unsigned int i = 0; i < order; ++i) {
+                    cg.insert(sub);
+                }
+                if (fac.mode == 1) cg.loop();
+                cg.pop();
+            }
+        }
+
+        return cg;
+    };
+
+    return parser;
+}
+
+#ifndef NDEBUG
+
+peg::parser build_ast_parser() {
+    peg::parser parser;
+    parser.set_logger([](size_t line, size_t col, const std::string &msg, const std::string &rule) {
+        fmt::print(stderr, "{}:{} [{}] {}\n", line, col, rule, msg);
+    });
+    auto ok = parser.load_grammar(GRAMMAR);
+    assert(ok);
+
+    parser.enable_ast();
+
+    return parser;
+}
+
+#endif
+
+std::vector<Op> compile(const std::string &source) {
+#ifndef NDEBUG
+    static peg::parser ast_parser = build_ast_parser();
+    std::shared_ptr<peg::Ast> ast;
+    bool ast_ok = ast_parser.parse(source, ast);
+    assert(ast_ok);
+//    std::cout << peg::ast_to_s(ast) << std::endl;
+#endif
+
+    static peg::parser parser = build_parser();
+    codegen cg;
+    bool ok = parser.parse(source, cg);
+    assert(ok);
+    return cg.ops;
+}
+
+Graph eval(const std::vector<Op> &ops) {
+    std::vector<std::stack<size_t>> stacks(1);
+
+    Graph g;
+    stacks.back().emplace(g.rank++);
+
+    for (const auto &op: ops) {
+        switch (op.code) {
+            case Op::FREE:
+            case Op::LINK: {
+                tc::Mult order = tc::FREE;
+
+                if (op.code == Op::LINK) {
+                    order = op.value;
+                }
+
+                auto top = stacks.back().top();
+                auto curr = g.rank++;
+
+                stacks.back().emplace(curr);
+                g.edges.emplace_back(top, curr, order);
+
+                break;
+            }
+            case Op::PUSH: {
+                stacks.emplace_back();
+
+                auto ptop = stacks[stacks.size() - 2].top();
+                stacks.back().emplace(ptop);
+
+                break;
+            }
+            case Op::POP: {
+                stacks.pop_back();
+
+                break;
+            }
+            case Op::LOOP: {
+                g.rank--;
+
+                auto ptop = stacks[stacks.size() - 2].top();
+                auto &[top, _, order] = g.edges.back();
+
+                stacks.back().pop();
+                g.edges.pop_back();
+
+                g.edges.emplace_back(top, ptop, order);
+
+                break;
+            }
+            default:
+                throw std::runtime_error("Invalid opcode");
+        }
+    }
+
+    return g;
+}
+
+namespace tc {
+    Group coxeter(const std::string &symbol) {
+        auto ops = compile(symbol);
+        auto diagram = eval(ops);
+        Group res(diagram.rank);
+        for (const auto &[i, j, m]: diagram.edges) {
+            res.set(i, j, m);
+        }
+        return res;
+    }
+}

tc/src/solve.cpp

@@ -0,0 +1,211 @@
+#include <algorithm>
+#include <queue>
+#include <utility>
+#include <vector>
+
+#include <tc/core.hpp>
+
+namespace tc {
+    /**
+     * Each coset is associated a row in each table.
+     * Rows document the "loops" formed by 
+     */
+    struct Row {
+        bool free: 1;
+        bool idem: 1;
+        unsigned int gnr: 14;  // progress through the loop
+        unsigned int lst_idx: 32;  // the coset that would complete the loop
+
+        Row() : free(true), idem(false), gnr(0), lst_idx(0) {}
+    };
+
+    struct Tables {
+        std::vector<Group::Rel> rels;
+        std::vector<std::vector<Row>> rows;
+
+        explicit Tables(std::vector<Group::Rel> rels)
+            : rels(std::move(rels)), rows() {
+        }
+
+        [[nodiscard]] size_t size() const {
+            return rels.size();
+        }
+
+        void add_row() {
+            rows.emplace_back(rels.size());
+        }
+    };
+
+    [[nodiscard]] Cosets Group::solve(std::vector<size_t> const &idxs, size_t bound) const {
+        // region Initialize Cosets Table
+        Cosets cosets(rank());
+        cosets.add_row();
+
+        if (rank() == 0) {
+            cosets._complete = true;
+            return cosets;
+        }
+
+        for (size_t g: idxs) {
+            if (g < rank())
+                cosets.set(0, g, 0);
+        }
+        // endregion
+
+        // region Initialize Relation Tables
+        std::vector<Group::Rel> rels;
+        for (int i = 0; i < rank(); ++i) {
+            for (int j = i + 1; j < rank(); ++j) {
+                // The algorithm only works for Coxeter groups; multiplicities m_ii=1 are assumed. Relation tables
+                // _may_ be added for them, but they are redundant and hurt performance so are skipped.
+                if (i == j) continue;
+
+                // Coxeter groups admit infinite multiplicities, represented by contexpr tc::FREE. Relation tables
+                // for these should be skipped.
+                auto m = get(i, j);
+
+                if (m == FREE) {
+                    continue;
+                }
+
+                rels.emplace_back(i, j, m);
+            }
+        }
+
+        Tables rel_tables(rels);
+        std::vector<std::vector<size_t>> tables_for(rank());
+        int rel_idx = 0;
+        for (const auto &[i, j, m]: rels) {
+            tables_for[i].push_back(rel_idx);
+            tables_for[j].push_back(rel_idx);
+            rel_idx++;
+        }
+
+        std::vector<size_t> lst_vals;
+        rel_tables.add_row();
+        for (int table_idx = 0; table_idx < rel_tables.size(); ++table_idx) {
+            const auto &[i, j, m] = rel_tables.rels[table_idx];
+            Row &row = rel_tables.rows[0][table_idx];
+
+            if (!cosets.isset(0, i) && !cosets.isset(0, j)) {
+                row.lst_idx = lst_vals.size();
+                lst_vals.push_back(0);
+                row.free = false;
+                row.gnr = 0;
+            } else {
+                row.free = false;
+                row.gnr = 1;
+                row.idem = true;
+            }
+        }
+        // endregion
+
+        size_t idx = 0;
+        size_t fact_idx;
+        size_t coset, gen, target, lst;
+
+        while (true) {
+            // find next unknown product
+            while (idx < cosets.size() and cosets.isset(idx))
+                idx++;
+
+            if (cosets.order() >= bound) {
+                return cosets;
+            }
+
+            // if there are none, then return
+            if (idx == cosets.size()) {
+                // todo unrolled linked list interval
+//                rel_tables.del_rows_to(idx / ngens);  
+                break;
+            }
+
+            // the unknown product must be a new coset, so add it
+            target = cosets.order();
+            cosets.add_row();
+            rel_tables.add_row();
+
+            // queue of products that equal target
+            std::queue<size_t> facts;
+            facts.push(idx);  // new product should be recorded and propagated
+
+            // todo unrolled linked list interval
+//            rel_tables.del_rows_to(coset);
+
+            // find all products which also lead to target
+            while (!facts.empty()) {
+                fact_idx = facts.front();
+                facts.pop();
+
+                // skip if this product was already learned
+                if (cosets.get(fact_idx) != -1) continue;
+
+                cosets.set(fact_idx, target);
+
+                coset = fact_idx / rank();
+                gen = fact_idx % rank();
+
+                // If the product stays within the coset todo
+                for (size_t table_idx: tables_for[gen]) {
+                    auto &[i, j, m] = rel_tables.rels[table_idx];
+                    auto &trow = rel_tables.rows[target][table_idx];
+                    auto &crow = rel_tables.rows[coset][table_idx];
+
+                    size_t other_gen = (i == gen) ? j : i;
+
+                    // Test if loop is closed
+                    if (trow.free) {
+                        trow = crow;
+                        trow.gnr++;
+
+                        if (target == coset) {
+                            trow.idem = true;
+                        }
+
+                        if (trow.idem) {
+                            if (trow.gnr == m) {
+                                // loop is closed, but idempotent, so the target links to itself via the other generator.
+                                // todo might be able to move this logic up into the (target == coset) block and avoid those computations.
+                                facts.push(target * rank() + other_gen);
+                            }
+                        } else {
+                            if (trow.gnr == m - 1) {
+                                // loop is almost closed. record that the target closes this loop.
+                                lst_vals[trow.lst_idx] = target;
+                            } else if (trow.gnr == m) {
+                                // loop is closed. We know the last element in the loop must link with this one. 
+                                lst = lst_vals[trow.lst_idx];
+//                            delete trow.lst_ptr;
+                                facts.push(lst * rank() + other_gen);
+                            }
+                        }
+                    }
+                }
+            }
+
+            // If any target row wasn't identified with a loop,
+            // then assign it a new loop.
+            for (size_t table_idx = 0; table_idx < rel_tables.size(); table_idx++) {
+                auto &[i, j, m] = rel_tables.rels[table_idx];
+                auto &trow = rel_tables.rows[target][table_idx];
+
+                if (trow.free) {
+                    if ((cosets.get(target, i) != target) and
+                        (cosets.get(target, j) != target)) {
+                        trow.lst_idx = lst_vals.size();
+                        trow.free = false;
+                        lst_vals.push_back(0);
+                        trow.gnr = 0;
+                    } else {
+                        trow.free = false;
+                        trow.gnr = 1;
+                        trow.idem = true;
+                    }
+                }
+            }
+        }
+
+        cosets._complete = true;
+        return cosets;
+    }
+}

tc/test/CMakeLists.txt

@@ -0,0 +1,18 @@
+include(GoogleTest)
+
+add_executable(test_solve test_solve.cpp)
+target_link_libraries(test_solve tc GTest::gtest_main)
+
+add_executable(test_lang test_lang.cpp)
+target_link_libraries(test_lang PUBLIC tc::tc GTest::gtest_main)
+
+set(MIN_DEBUG_CPS 200000)
+set(MIN_RELEASE_CPS 1000000)
+
+target_compile_definitions(
+    test_solve PUBLIC
+    MIN_COS_PER_SEC=$<IF:$<CONFIG:Debug>,${MIN_DEBUG_CPS},${MIN_RELEASE_CPS}>
+)
+
+gtest_discover_tests(test_solve)
+gtest_discover_tests(test_lang)

tc/test/test_lang.cpp

@@ -0,0 +1,54 @@
+#include <tc/core.hpp>
+#include <tc/groups.hpp>
+
+#include <gtest/gtest.h>
+
+// "5 3 3"
+// "5 (3 3)"
+// "[5 3 3]"
+// "[4 3 [3 5] 3]"
+// "{3 4 5 6 7 8 9}"
+// "3 {3 3 [4] 3} 5"
+// "5 * 3"
+// "5 * [3]"
+// "5 * {2 3}"
+// "5 * [3 2]"
+// "(5 2) * [3 2]"
+// "4 [3 * [2 3]] 5"
+// "{3 * 3} [4] [5]"
+
+TEST(coxeter, simple) {
+    auto g = tc::coxeter("5 3 3");
+
+    ASSERT_EQ(g.rank(), 4);
+
+    EXPECT_EQ(g.get(0, 0), 1);
+    EXPECT_EQ(g.get(0, 1), 5);
+    EXPECT_EQ(g.get(0, 2), 2);
+    EXPECT_EQ(g.get(0, 3), 2);
+    
+    EXPECT_EQ(g.get(1, 0), 5);
+    EXPECT_EQ(g.get(1, 1), 1);
+    EXPECT_EQ(g.get(1, 2), 3);
+    EXPECT_EQ(g.get(1, 3), 2);
+    
+    EXPECT_EQ(g.get(2, 0), 2);
+    EXPECT_EQ(g.get(2, 1), 3);
+    EXPECT_EQ(g.get(2, 2), 1);
+    EXPECT_EQ(g.get(2, 3), 3);
+    
+    EXPECT_EQ(g.get(3, 0), 2);
+    EXPECT_EQ(g.get(3, 1), 2);
+    EXPECT_EQ(g.get(3, 2), 3);
+    EXPECT_EQ(g.get(3, 3), 1);
+}
+
+TEST(coxeter, looping) {
+    auto g = tc::coxeter("{5 3 4}");
+
+    ASSERT_EQ(g.rank(), 3);
+
+    EXPECT_EQ(g.get(0, 1), 5);
+    EXPECT_EQ(g.get(1, 2), 3);
+    EXPECT_EQ(g.get(2, 0), 4);
+}

tc/test/test_solve.cpp

@@ -0,0 +1,225 @@
+#include <ctime>
+#include <vector>
+
+#include <tc/groups.hpp>
+#include <tc/core.hpp>
+
+#include <gtest/gtest.h>
+
+/// helper for testing solve and testing speed
+testing::AssertionResult AssertSolveOrder(
+    const char *group_expr,
+    const char *sub_gens_expr,
+    const char *expected_order_expr,
+    const tc::Group &group,
+    const std::vector<size_t> &sub_gens,
+    size_t expected_order
+) {
+    auto start = std::clock();
+    auto cosets = group.solve(sub_gens);
+    auto end = std::clock();
+
+    size_t actual_order = cosets.order();
+
+    auto total_sec = (double) (end - start) / CLOCKS_PER_SEC;
+    auto cosets_per_sec = (double) actual_order / total_sec;
+
+    bool order_good = actual_order == expected_order;
+    bool speed_good = cosets_per_sec >= MIN_COS_PER_SEC || total_sec < 0.0001;
+    // extremely short times cause false negatives. ex. A2 can be solved in only 3 clocks.
+
+    if (order_good && speed_good) {
+        return testing::AssertionSuccess();
+    }
+
+    testing::AssertionResult res = testing::AssertionFailure();
+    res << group_expr << " / " << sub_gens_expr << " :";
+    if (!order_good) {
+        res << " Gave order " << actual_order << " but expected order " << expected_order << ".";
+    }
+    if (!speed_good) {
+        res << " Solution too slow (" << cosets_per_sec
+            << " cos/s < " << MIN_COS_PER_SEC << ")."
+            << " " << std::fixed << total_sec << " s.";
+    }
+
+    return res;
+}
+
+#define EXPECT_SOLVE_ORDER(group, sub_gens, expected_order) EXPECT_PRED_FORMAT3(AssertSolveOrder, group, sub_gens, expected_order);
+
+using v = std::vector<size_t>;
+
+tc::Group A(unsigned int n) {
+    return tc::vcoxeter("3 * {}", {n - 1});
+}
+
+tc::Group B(unsigned int n) {
+    return tc::vcoxeter("4 3 * {}", {n - 2});
+}
+
+tc::Group D(unsigned int n) {
+    return tc::vcoxeter("3 * [1 1 {}]", {n - 3});
+}
+
+tc::Group E(unsigned int n) {
+    return tc::vcoxeter("3 * [1 2 {}]", {n - 4});
+}
+
+tc::Group F4() {
+    return tc::coxeter("3 4 3");
+}
+
+tc::Group G2() {
+    return tc::coxeter("6");
+}
+
+tc::Group H(unsigned int n) {
+    return tc::vcoxeter("5 3 * {}", {n - 2});
+}
+
+tc::Group I2(unsigned int n) {
+    return tc::vcoxeter("{}", {n});
+}
+
+tc::Group T(unsigned int m, unsigned int n) {
+    return tc::vcoxeter("{} 2 {}", {m, n});
+}
+
+tc::Group T(unsigned int n) {
+    return T(n, n);
+}
+
+// See the group orders here https://en.wikipedia.org/wiki/Coxeter_group#Properties
+
+TEST(solve, A) {
+    EXPECT_SOLVE_ORDER(A(1), v({}), 2);
+    EXPECT_SOLVE_ORDER(A(2), v({}), 6);
+    EXPECT_SOLVE_ORDER(A(3), v({}), 24);
+    EXPECT_SOLVE_ORDER(A(3), v({0}), 12);
+    EXPECT_SOLVE_ORDER(A(3), v({0, 1}), 4);
+    EXPECT_SOLVE_ORDER(A(3), v({0, 2}), 6);
+    EXPECT_SOLVE_ORDER(A(3), v({2}), 12);
+    EXPECT_SOLVE_ORDER(A(4), v({}), 120);
+    EXPECT_SOLVE_ORDER(A(4), v({0}), 60);
+    EXPECT_SOLVE_ORDER(A(4), v({0, 1}), 20);
+    EXPECT_SOLVE_ORDER(A(4), v({2}), 60);
+    EXPECT_SOLVE_ORDER(A(4), v({0, 2}), 30);
+}
+
+TEST(solve, B) {
+    EXPECT_SOLVE_ORDER(B(2), v({}), 8);
+    EXPECT_SOLVE_ORDER(B(3), v({}), 48);
+    EXPECT_SOLVE_ORDER(B(3), v({0}), 24);
+    EXPECT_SOLVE_ORDER(B(3), v({0, 2}), 12);
+    EXPECT_SOLVE_ORDER(B(4), v({}), 384);
+    EXPECT_SOLVE_ORDER(B(4), v({0}), 192);
+    EXPECT_SOLVE_ORDER(B(4), v({0, 2}), 96);
+    EXPECT_SOLVE_ORDER(B(5), v({}), 3840);
+    EXPECT_SOLVE_ORDER(B(5), v({0}), 1920);
+    EXPECT_SOLVE_ORDER(B(5), v({0, 2}), 960);
+    EXPECT_SOLVE_ORDER(B(5), v({0, 2, 3}), 320);
+    EXPECT_SOLVE_ORDER(B(6), v({}), 46080);
+    EXPECT_SOLVE_ORDER(B(6), v({0}), 23040);
+    EXPECT_SOLVE_ORDER(B(6), v({0, 2}), 11520);
+    EXPECT_SOLVE_ORDER(B(6), v({0, 2, 3}), 3840);
+    EXPECT_SOLVE_ORDER(B(6), v({0, 2, 3, 5}), 1920);
+}
+
+TEST(solve, D) {
+    EXPECT_SOLVE_ORDER(D(3), v({}), 24);
+    EXPECT_SOLVE_ORDER(D(4), v({}), 192);
+    EXPECT_SOLVE_ORDER(D(4), v({0, 1}), 32);
+    EXPECT_SOLVE_ORDER(D(4), v({0, 1, 3}), 8);
+    EXPECT_SOLVE_ORDER(D(5), v({}), 1920);
+    EXPECT_SOLVE_ORDER(D(5), v({0, 1}), 320);
+    EXPECT_SOLVE_ORDER(D(5), v({0, 1, 3}), 80);
+    EXPECT_SOLVE_ORDER(D(5), v({0, 1, 3, 4}), 16);
+    EXPECT_SOLVE_ORDER(D(6), v({}), 23040);
+    EXPECT_SOLVE_ORDER(D(6), v({0, 1}), 3840);
+    EXPECT_SOLVE_ORDER(D(6), v({0, 1, 3}), 960);
+    EXPECT_SOLVE_ORDER(D(6), v({0, 1, 3, 5}), 480);
+}
+
+TEST(solve, E) {
+    EXPECT_SOLVE_ORDER(E(4), v({}), 120);
+    EXPECT_SOLVE_ORDER(E(4), v({2}), 60);
+    EXPECT_SOLVE_ORDER(E(4), v({2, 1}), 30);
+    EXPECT_SOLVE_ORDER(E(4), v({2, 1, 3}), 10);
+    EXPECT_SOLVE_ORDER(E(5), v({}), 1920);
+    EXPECT_SOLVE_ORDER(E(5), v({2}), 960);
+    EXPECT_SOLVE_ORDER(E(5), v({2, 1}), 480);
+    EXPECT_SOLVE_ORDER(E(5), v({2, 1, 3}), 160);
+    EXPECT_SOLVE_ORDER(E(6), v({}), 51840);
+    EXPECT_SOLVE_ORDER(E(6), v({2}), 25920);
+    EXPECT_SOLVE_ORDER(E(6), v({2, 1}), 12960);
+    EXPECT_SOLVE_ORDER(E(6), v({2, 1, 3}), 4320);
+}
+
+TEST(solve, F) {
+    EXPECT_SOLVE_ORDER(F4(), v({}), 1152);
+    EXPECT_SOLVE_ORDER(F4(), v({0}), 576);
+    EXPECT_SOLVE_ORDER(F4(), v({0, 2}), 288);
+    EXPECT_SOLVE_ORDER(F4(), v({1, 3}), 288);
+    EXPECT_SOLVE_ORDER(F4(), v({1, 2, 3}), 24);
+}
+
+TEST(solve, G) {
+    EXPECT_SOLVE_ORDER(G2(), v({}), 12);
+    EXPECT_SOLVE_ORDER(G2(), v({0}), 6);
+    EXPECT_SOLVE_ORDER(G2(), v({1}), 6);
+}
+
+TEST(solve, H) {
+    EXPECT_SOLVE_ORDER(H(2), v({}), 10);
+    EXPECT_SOLVE_ORDER(H(2), v({0}), 5);
+    EXPECT_SOLVE_ORDER(H(2), v({1}), 5);
+    EXPECT_SOLVE_ORDER(H(3), v({}), 120);
+    EXPECT_SOLVE_ORDER(H(3), v({0}), 60);
+    EXPECT_SOLVE_ORDER(H(3), v({0, 1}), 12);
+    EXPECT_SOLVE_ORDER(H(3), v({0, 2}), 30);
+    EXPECT_SOLVE_ORDER(H(3), v({1, 2}), 20);
+    EXPECT_SOLVE_ORDER(H(4), v({}), 14400);
+    EXPECT_SOLVE_ORDER(H(4), v({0}), 7200);
+    EXPECT_SOLVE_ORDER(H(4), v({1}), 7200);
+    EXPECT_SOLVE_ORDER(H(4), v({1, 3}), 3600);
+}
+
+TEST(solve, I) {
+    EXPECT_SOLVE_ORDER(I2(2), v({}), 4);
+    EXPECT_SOLVE_ORDER(I2(3), v({}), 6);
+    EXPECT_SOLVE_ORDER(I2(3), v({0}), 3);
+    EXPECT_SOLVE_ORDER(I2(3), v({1}), 3);
+    EXPECT_SOLVE_ORDER(I2(4), v({}), 8);
+    EXPECT_SOLVE_ORDER(I2(4), v({0}), 4);
+    EXPECT_SOLVE_ORDER(I2(4), v({1}), 4);
+    EXPECT_SOLVE_ORDER(I2(5), v({}), 10);
+    EXPECT_SOLVE_ORDER(I2(5), v({0}), 5);
+    EXPECT_SOLVE_ORDER(I2(5), v({1}), 5);
+}
+
+TEST(solve, T) {
+    EXPECT_SOLVE_ORDER(T(3), v({}), 36);
+    EXPECT_SOLVE_ORDER(T(4), v({}), 64);
+    EXPECT_SOLVE_ORDER(T(400), v({}), 640000);
+    EXPECT_SOLVE_ORDER(T(400), v({0}), 320000);
+    EXPECT_SOLVE_ORDER(T(400), v({0, 2}), 160000);
+    EXPECT_SOLVE_ORDER(T(400, 300), v({}), 480000);
+    EXPECT_SOLVE_ORDER(T(400, 300), v({0}), 240000);
+    EXPECT_SOLVE_ORDER(T(400, 300), v({0, 2}), 120000);
+}
+
+TEST(solve_large, B) {
+    EXPECT_SOLVE_ORDER(B(7), v({}), 645120);
+    EXPECT_SOLVE_ORDER(B(8), v({}), 10321920);
+}
+
+TEST(solve_large, E) {
+    EXPECT_SOLVE_ORDER(E(6), v({}), 51840);
+    EXPECT_SOLVE_ORDER(E(7), v({}), 2903040);
+}
+
+TEST(solve_large, T) {
+    EXPECT_SOLVE_ORDER(T(500), v({}), 1000000);
+    EXPECT_SOLVE_ORDER(T(1000), v({}), 4000000);
+}

vendor/glad/CMakeLists.txt

@@ -1,2 +0,0 @@
-add_library(glad src/glad.c)
-target_include_directories(glad PUBLIC include)

vendor/glad/include/KHR/khrplatform.h

@@ -1,290 +0,0 @@
-#ifndef __khrplatform_h_
-#define __khrplatform_h_
-
-/*
-** Copyright (c) 2008-2018 The Khronos Group Inc.
-**
-** Permission is hereby granted, free of charge, to any person obtaining a
-** copy of this software and/or associated documentation files (the
-** "Materials"), to deal in the Materials without restriction, including
-** without limitation the rights to use, copy, modify, merge, publish,
-** distribute, sublicense, and/or sell copies of the Materials, and to
-** permit persons to whom the Materials are furnished to do so, subject to
-** the following conditions:
-**
-** The above copyright notice and this permission notice shall be included
-** in all copies or substantial portions of the Materials.
-**
-** THE MATERIALS ARE PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
-** EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF
-** MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT.
-** IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY
-** CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT,
-** TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE
-** MATERIALS OR THE USE OR OTHER DEALINGS IN THE MATERIALS.
-*/
-
-/* Khronos platform-specific types and definitions.
- *
- * The master copy of khrplatform.h is maintained in the Khronos EGL
- * Registry repository at https://github.com/KhronosGroup/EGL-Registry
- * The last semantic modification to khrplatform.h was at commit ID:
- *      67a3e0864c2d75ea5287b9f3d2eb74a745936692
- *
- * Adopters may modify this file to suit their platform. Adopters are
- * encouraged to submit platform specific modifications to the Khronos
- * group so that they can be included in future versions of this file.
- * Please submit changes by filing pull requests or issues on
- * the EGL Registry repository linked above.
- *
- *
- * See the Implementer's Guidelines for information about where this file
- * should be located on your system and for more details of its use:
- *    http://www.khronos.org/registry/implementers_guide.pdf
- *
- * This file should be included as
- *        #include <KHR/khrplatform.h>
- * by Khronos client API header files that use its types and defines.
- *
- * The types in khrplatform.h should only be used to define API-specific types.
- *
- * Types defined in khrplatform.h:
- *    khronos_int8_t              signed   8  bit
- *    khronos_uint8_t             unsigned 8  bit
- *    khronos_int16_t             signed   16 bit
- *    khronos_uint16_t            unsigned 16 bit
- *    khronos_int32_t             signed   32 bit
- *    khronos_uint32_t            unsigned 32 bit
- *    khronos_int64_t             signed   64 bit
- *    khronos_uint64_t            unsigned 64 bit
- *    khronos_intptr_t            signed   same number of bits as a pointer
- *    khronos_uintptr_t           unsigned same number of bits as a pointer
- *    khronos_ssize_t             signed   size
- *    khronos_usize_t             unsigned size
- *    khronos_float_t             signed   32 bit floating point
- *    khronos_time_ns_t           unsigned 64 bit time in nanoseconds
- *    khronos_utime_nanoseconds_t unsigned time interval or absolute time in
- *                                         nanoseconds
- *    khronos_stime_nanoseconds_t signed time interval in nanoseconds
- *    khronos_boolean_enum_t      enumerated boolean type. This should
- *      only be used as a base type when a client API's boolean type is
- *      an enum. Client APIs which use an integer or other type for
- *      booleans cannot use this as the base type for their boolean.
- *
- * Tokens defined in khrplatform.h:
- *
- *    KHRONOS_FALSE, KHRONOS_TRUE Enumerated boolean false/true values.
- *
- *    KHRONOS_SUPPORT_INT64 is 1 if 64 bit integers are supported; otherwise 0.
- *    KHRONOS_SUPPORT_FLOAT is 1 if floats are supported; otherwise 0.
- *
- * Calling convention macros defined in this file:
- *    KHRONOS_APICALL
- *    KHRONOS_APIENTRY
- *    KHRONOS_APIATTRIBUTES
- *
- * These may be used in function prototypes as:
- *
- *      KHRONOS_APICALL void KHRONOS_APIENTRY funcname(
- *                                  int arg1,
- *                                  int arg2) KHRONOS_APIATTRIBUTES;
- */
-
-#if defined(__SCITECH_SNAP__) && !defined(KHRONOS_STATIC)
-#   define KHRONOS_STATIC 1
-#endif
-
-/*-------------------------------------------------------------------------
- * Definition of KHRONOS_APICALL
- *-------------------------------------------------------------------------
- * This precedes the return type of the function in the function prototype.
- */
-#if defined(KHRONOS_STATIC)
-    /* If the preprocessor constant KHRONOS_STATIC is defined, make the
-     * header compatible with static linking. */
-#   define KHRONOS_APICALL
-#elif defined(_WIN32)
-#   define KHRONOS_APICALL __declspec(dllimport)
-#elif defined (__SYMBIAN32__)
-#   define KHRONOS_APICALL IMPORT_C
-#elif defined(__ANDROID__)
-#   define KHRONOS_APICALL __attribute__((visibility("default")))
-#else
-#   define KHRONOS_APICALL
-#endif
-
-/*-------------------------------------------------------------------------
- * Definition of KHRONOS_APIENTRY
- *-------------------------------------------------------------------------
- * This follows the return type of the function  and precedes the function
- * name in the function prototype.
- */
-#if defined(_WIN32) && !defined(_WIN32_WCE) && !defined(KHRONOS_STATIC)
-    /* Win32 but not WinCE */
-#   define KHRONOS_APIENTRY __stdcall
-#else
-#   define KHRONOS_APIENTRY
-#endif
-
-/*-------------------------------------------------------------------------
- * Definition of KHRONOS_APIATTRIBUTES
- *-------------------------------------------------------------------------
- * This follows the closing parenthesis of the function prototype arguments.
- */
-#if defined (__ARMCC_2__)
-#define KHRONOS_APIATTRIBUTES __softfp
-#else
-#define KHRONOS_APIATTRIBUTES
-#endif
-
-/*-------------------------------------------------------------------------
- * basic type definitions
- *-----------------------------------------------------------------------*/
-#if (defined(__STDC_VERSION__) && __STDC_VERSION__ >= 199901L) || defined(__GNUC__) || defined(__SCO__) || defined(__USLC__)
-
-
-/*
- * Using <stdint.h>
- */
-#include <stdint.h>
-typedef int32_t                 khronos_int32_t;
-typedef uint32_t                khronos_uint32_t;
-typedef int64_t                 khronos_int64_t;
-typedef uint64_t                khronos_uint64_t;
-#define KHRONOS_SUPPORT_INT64   1
-#define KHRONOS_SUPPORT_FLOAT   1
-
-#elif defined(__VMS ) || defined(__sgi)
-
-/*
- * Using <inttypes.h>
- */
-#include <inttypes.h>
-typedef int32_t                 khronos_int32_t;
-typedef uint32_t                khronos_uint32_t;
-typedef int64_t                 khronos_int64_t;
-typedef uint64_t                khronos_uint64_t;
-#define KHRONOS_SUPPORT_INT64   1
-#define KHRONOS_SUPPORT_FLOAT   1
-
-#elif defined(_WIN32) && !defined(__SCITECH_SNAP__)
-
-/*
- * Win32
- */
-typedef __int32                 khronos_int32_t;
-typedef unsigned __int32        khronos_uint32_t;
-typedef __int64                 khronos_int64_t;
-typedef unsigned __int64        khronos_uint64_t;
-#define KHRONOS_SUPPORT_INT64   1
-#define KHRONOS_SUPPORT_FLOAT   1
-
-#elif defined(__sun__) || defined(__digital__)
-
-/*
- * Sun or Digital
- */
-typedef int                     khronos_int32_t;
-typedef unsigned int            khronos_uint32_t;
-#if defined(__arch64__) || defined(_LP64)
-typedef long int                khronos_int64_t;
-typedef unsigned long int       khronos_uint64_t;
-#else
-typedef long long int           khronos_int64_t;
-typedef unsigned long long int  khronos_uint64_t;
-#endif /* __arch64__ */
-#define KHRONOS_SUPPORT_INT64   1
-#define KHRONOS_SUPPORT_FLOAT   1
-
-#elif 0
-
-/*
- * Hypothetical platform with no float or int64 support
- */
-typedef int                     khronos_int32_t;
-typedef unsigned int            khronos_uint32_t;
-#define KHRONOS_SUPPORT_INT64   0
-#define KHRONOS_SUPPORT_FLOAT   0
-
-#else
-
-/*
- * Generic fallback
- */
-#include <stdint.h>
-typedef int32_t                 khronos_int32_t;
-typedef uint32_t                khronos_uint32_t;
-typedef int64_t                 khronos_int64_t;
-typedef uint64_t                khronos_uint64_t;
-#define KHRONOS_SUPPORT_INT64   1
-#define KHRONOS_SUPPORT_FLOAT   1
-
-#endif
-
-
-/*
- * Types that are (so far) the same on all platforms
- */
-typedef signed   char          khronos_int8_t;
-typedef unsigned char          khronos_uint8_t;
-typedef signed   short int     khronos_int16_t;
-typedef unsigned short int     khronos_uint16_t;
-
-/*
- * Types that differ between LLP64 and LP64 architectures - in LLP64,
- * pointers are 64 bits, but 'long' is still 32 bits. Win64 appears
- * to be the only LLP64 architecture in current use.
- */
-#ifdef _WIN64
-typedef signed   long long int khronos_intptr_t;
-typedef unsigned long long int khronos_uintptr_t;
-typedef signed   long long int khronos_ssize_t;
-typedef unsigned long long int khronos_usize_t;
-#else
-typedef signed   long  int     khronos_intptr_t;
-typedef unsigned long  int     khronos_uintptr_t;
-typedef signed   long  int     khronos_ssize_t;
-typedef unsigned long  int     khronos_usize_t;
-#endif
-
-#if KHRONOS_SUPPORT_FLOAT
-/*
- * Float type
- */
-typedef          float         khronos_float_t;
-#endif
-
-#if KHRONOS_SUPPORT_INT64
-/* Time types
- *
- * These types can be used to represent a time interval in nanoseconds or
- * an absolute Unadjusted System Time.  Unadjusted System Time is the number
- * of nanoseconds since some arbitrary system event (e.g. since the last
- * time the system booted).  The Unadjusted System Time is an unsigned
- * 64 bit value that wraps back to 0 every 584 years.  Time intervals
- * may be either signed or unsigned.
- */
-typedef khronos_uint64_t       khronos_utime_nanoseconds_t;
-typedef khronos_int64_t        khronos_stime_nanoseconds_t;
-#endif
-
-/*
- * Dummy value used to pad enum types to 32 bits.
- */
-#ifndef KHRONOS_MAX_ENUM
-#define KHRONOS_MAX_ENUM 0x7FFFFFFF
-#endif
-
-/*
- * Enumerated boolean type
- *
- * Values other than zero should be considered to be true.  Therefore
- * comparisons should not be made against KHRONOS_TRUE.
- */
-typedef enum {
-    KHRONOS_FALSE = 0,
-    KHRONOS_TRUE  = 1,
-    KHRONOS_BOOLEAN_ENUM_FORCE_SIZE = KHRONOS_MAX_ENUM
-} khronos_boolean_enum_t;
-
-#endif /* __khrplatform_h_ */

vis/CMakeLists.txt

@@ -1,9 +1,4 @@
-add_custom_target(shaders ALL DEPENDS shader_output)
-add_custom_command(
-        OUTPUT shader_output
-        COMMAND ${CMAKE_COMMAND} -E copy_directory ${CMAKE_CURRENT_SOURCE_DIR}/shaders ${CMAKE_CURRENT_BINARY_DIR}/shaders
-        COMMENT "copied shaders"
-)
+add_subdirectory(shaders)
 
 add_custom_target(presets ALL DEPENDS preset_output)
 add_custom_command(
@@ -12,10 +7,13 @@ add_custom_command(
         COMMENT "copied preses"
 )
 
-add_library(vis-util INTERFACE)
-target_include_directories(vis-util INTERFACE include)
-
-add_executable(vis src/main.cpp)
+add_executable(vis
+    src/main.cpp
+    src/components.hpp
+    )
 target_include_directories(vis PRIVATE include)
-target_link_libraries(vis PRIVATE tc glad glm glfw yaml-cpp)
-add_dependencies(vis shaders presets)
+target_link_libraries(vis PRIVATE
+    tc shaders
+    glad eigen glfw yaml-cpp fmt EnTT imgui
+    )
+add_dependencies(vis presets)

vis/include/cgl/buffer.hpp

@@ -8,22 +8,15 @@ namespace cgl {
     template<class T>
     class Buffer {
         GLuint id{};
+        size_t _count = 0;
 
     public:
+        using Element = T;
+
         Buffer() {
             glCreateBuffers(1, &id);
         }
 
-        Buffer(const T &data, GLenum usage = GL_STATIC_DRAW)
-            : Buffer() {
-            put(data, usage);
-        }
-
-        Buffer(const std::vector<T> &data, GLenum usage = GL_STATIC_DRAW)
-            : Buffer() {
-            put(data, usage);
-        }
-
         Buffer(Buffer &) = delete;
 
         Buffer(Buffer &&o) noexcept {
@@ -40,28 +33,26 @@ namespace cgl {
         }
 
         [[nodiscard]] size_t size() const {
-            GLint res;
-            glGetNamedBufferParameteriv(id, GL_BUFFER_SIZE, &res);
-            return (size_t) res;
+            return _count * sizeof(T);
         }
 
         [[nodiscard]] size_t count() const {
-            return size() / sizeof(T);
+            return _count;
+        }
+
+        template<typename It>
+        void put(It begin, It end, GLenum usage = GL_STATIC_DRAW) {
+            _count = end - begin;
+
+            glNamedBufferData(id, sizeof(T) * _count, nullptr, usage);
+            void* ptr = glMapNamedBuffer(id, GL_WRITE_ONLY);
+            std::copy(begin, end, (T*) ptr);
+            glUnmapNamedBuffer(id);
         }
 
         void put(const T &data, GLenum usage = GL_STATIC_DRAW) {
-            glNamedBufferData(id, sizeof(T), &data, usage);
-        }
-
-        void put(const std::vector<T> &data, GLenum usage = GL_STATIC_DRAW) {
-            glNamedBufferData(id, sizeof(T) * data.size(), &data[0], usage);
-        }
-
-        void bound(GLenum target, const std::function<void()> &action) const {
-            glBindBuffer(target, id);
-            action();
-            glBindBuffer(target, 0);
+            T const* ptr = &data;
+            put(ptr, ptr + 1, usage);
         }
     };
-
 }

vis/include/cgl/debug.hpp

@@ -0,0 +1,85 @@
+#pragma once
+
+#include <glad/glad.h>
+#include <fmt/core.h>
+
+#include <string>
+
+#ifndef NDEBUG
+
+void GLAPIENTRY log_gl_debug_callback(
+    GLenum source,
+    GLenum type,
+    GLuint id,
+    GLenum severity,
+    GLsizei length,
+    const GLchar *message,
+    const void *userParam
+) {
+    std::string s_source;
+    switch (type) {
+        case GL_DEBUG_SOURCE_API:
+            s_source = "API:";
+        case GL_DEBUG_SOURCE_WINDOW_SYSTEM:
+            s_source = "WINDOW:";
+        case GL_DEBUG_SOURCE_SHADER_COMPILER:
+            s_source = "SHADER:";
+        case GL_DEBUG_SOURCE_THIRD_PARTY:
+            s_source = "3P:";
+        case GL_DEBUG_SOURCE_APPLICATION:
+            s_source = "APP:";
+        default:
+            s_source = "";
+    }
+
+    std::string s_type;
+    switch (type) {
+        case GL_DEBUG_TYPE_ERROR:
+            s_type = "ERROR:";
+            break;
+        case GL_DEBUG_TYPE_DEPRECATED_BEHAVIOR:
+            s_type = "DEPRECATED:";
+            break;
+        case GL_DEBUG_TYPE_UNDEFINED_BEHAVIOR:
+            s_type = "UNDEFINED:";
+            break;
+        case GL_DEBUG_TYPE_PORTABILITY:
+            s_type = "PORTABILITY:";
+            break;
+        case GL_DEBUG_TYPE_PERFORMANCE:
+            s_type = "PERFORMANCE:";
+            break;
+        case GL_DEBUG_TYPE_MARKER:
+            s_type = "MARKER:";
+            break;
+        case GL_DEBUG_TYPE_PUSH_GROUP:
+            s_type = "PUSH_GROUP:";
+            break;
+        case GL_DEBUG_TYPE_POP_GROUP:
+            s_type = "POP_GROUP:";
+            break;
+        default:
+            s_type = "";
+            break;
+    }
+
+    std::string s_severity;
+    switch (severity) {
+        case GL_DEBUG_SEVERITY_HIGH:
+            s_severity = "HIGH:";
+            break;
+        case GL_DEBUG_SEVERITY_MEDIUM:
+            s_severity = "MED:";
+            break;
+        case GL_DEBUG_SEVERITY_LOW:
+            s_severity = "LOW:";
+            break;
+        default:
+            s_severity = "INFO:";
+            break;
+    }
+
+    fmt::print(stderr, "GL:{}{}{} {}\n", s_source, s_type, s_severity, message);
+}
+
+#endif

vis/include/cgl/pipeline.hpp

@@ -83,11 +83,5 @@ namespace cgl {
             glUseProgramStages(id, GL_COMPUTE_SHADER_BIT, pgm);
             return *this;
         }
-
-        void bound(const std::function<void()> &action) const {
-            glBindProgramPipeline(id);
-            action();
-            glBindProgramPipeline(0);
-        }
     };
 }

vis/include/cgl/shaderprogram.hpp

@@ -19,7 +19,9 @@ namespace cgl{
             glProgramParameteri(id, GL_PROGRAM_SEPARABLE, GL_TRUE);
         }
 
-        ShaderProgram(const std::string &src) : ShaderProgram() {
+        explicit ShaderProgram(const std::string_view &src) : ShaderProgram(std::string(src)) {}
+
+        explicit ShaderProgram(const std::string &src) : ShaderProgram() {
             Shader<mode> sh(src);
 
             attach(sh);

vis/include/cgl/vertexarray.hpp

@@ -33,43 +33,34 @@ namespace cgl {
             return id;
         }
 
-        void bound(const std::function<void()> &action) const {
-            glBindVertexArray(id);
-            action();
-            glBindVertexArray(0);
-        }
-
-        template<class T>
-        void pointer(
+        void format(
             GLuint index,
-            const Buffer<T> &buf,
             unsigned size,
             GLenum type,
             bool normalized = false,
             unsigned stride = 0
-        ) const {
-            bound([&]() {
-                glEnableVertexAttribArray(index);
-                buf.bound(GL_ARRAY_BUFFER, [&]() {
-                    glVertexAttribPointer(index, size, type, normalized, stride, nullptr);
-                });
-            });
+        ) {
+            glEnableVertexArrayAttrib(id, index);
+            glVertexArrayAttribFormat(id, index, size, type, normalized, stride);
         }
 
-        template<class T>
-        void ipointer(
+        void iformat(
             GLuint index,
-            const Buffer<T> &buf,
             unsigned size,
             GLenum type,
             unsigned stride = 0
-        ) const {
-            bound([&]() {
-                glEnableVertexAttribArray(index);
-                buf.bound(GL_ARRAY_BUFFER, [&]() {
-                    glVertexAttribIPointer(index, size, type, stride, nullptr);
-                });
-            });
+        ) {
+            glEnableVertexArrayAttrib(id, index);
+            glVertexArrayAttribIFormat(id, index, size, type, stride);
+        }
+
+        template<class Buf>
+        void vertexBuffer(
+            GLuint index,
+            Buf &buf,
+            unsigned offset = 0
+        ) {
+            glVertexArrayVertexBuffer(id, index, buf, offset, sizeof(typename Buf::Element));
         }
     };
 }

vis/include/combo_iterator.hpp

@@ -1,98 +0,0 @@
-#pragma once
-
-#include <vector>
-#include <array>
-#include <algorithm>
-#include <numeric>
-#include <stdexcept>
-
-size_t choose(size_t n, size_t k) {
-    if (k == 0) return 1;
-    return n * choose(n - 1, k - 1) / k;
-}
-
-template<class T>
-class ComboIterator {
-private:
-    const std::vector<T> &options;
-
-    std::vector<bool> bits;
-    std::vector<T> curr;
-    int at;
-
-    void set_curr() {
-        for (int i = 0, j = 0; i < bits.size(); ++i) {
-            if (bits[i]) curr[j++] = options[i];
-        }
-    }
-
-public:
-    ComboIterator(const std::vector<T> &options, int k, int at = 0)
-        : options(options), bits(options.size()), curr(k), at(at) {
-        std::fill(bits.begin(), bits.begin() + k, true);
-        set_curr();
-    }
-
-    [[nodiscard]] bool operator==(const ComboIterator<T> &o) const {
-        return at == o.at;
-    }
-
-    [[nodiscard]] bool operator!=(const ComboIterator<T> &o) const {
-        return at != o.at;
-    }
-
-    auto operator*() const {
-        return curr;
-    }
-
-    const auto &operator->() const {
-        return &this;
-    }
-
-    auto operator++(int) {
-        std::prev_permutation(bits.begin(), bits.end());
-        set_curr();
-        ++at;
-        return *this;
-    }
-
-    auto operator++() &{
-        auto res = *this;
-        (*this)++;
-        return res;
-    }
-
-    auto operator--(int) {
-        std::next_permutation(bits.begin(), bits.end());
-        set_curr();
-        --at;
-        return *this;
-    }
-
-    auto operator--() &{
-        auto res = *this;
-        (*this)--;
-        return res;
-    }
-};
-
-template<class T>
-class Combos {
-private:
-    const std::vector<T> options;
-    int k;
-    int size;
-
-public:
-    Combos(const std::vector<T> &options, int k)
-        : options(options), k(k), size(choose(options.size(), k)) {
-    }
-
-    ComboIterator<T> begin() const {
-        return ComboIterator<T>(options, k);
-    }
-
-    ComboIterator<T> end() const {
-        return ComboIterator<T>(options, k, size);
-    }
-};

vis/include/geometry.hpp

@@ -5,292 +5,57 @@
 #include <optional>
 #include <numeric>
 #include <iostream>
-#include "combo_iterator.hpp"
 
-/**
- * An primitive stage N indices.
- * @tparam N
- */
+#include <Eigen/Eigen>
+
 template<unsigned N>
-struct Primitive {
-    static_assert(N > 0, "Primitives must contain at least one point. Primitive<0> or lower is impossible.");
+using Indices = Eigen::Array<unsigned, N, Eigen::Dynamic>;
 
-    std::array<unsigned, N> inds;
+template<int N>
+using vec = Eigen::Matrix<float, N, 1>;
+template<int N>
+using mat = Eigen::Matrix<float, N, N>;
 
-    Primitive() = default;
+using vec1 = vec<1>;
+using vec2 = vec<2>;
+using vec3 = vec<3>;
+using vec4 = vec<4>;
+using vec5 = vec<5>;
 
-    Primitive(const Primitive<N> &) = default;
-
-    Primitive(const Primitive<N - 1> &sub, unsigned root) {
-        std::copy(sub.inds.begin(), sub.inds.end(), inds.begin());
-        inds[N - 1] = root;
-    }
-
-    ~Primitive() = default;
-
-    inline void flip() {
-        if (N > 1) std::swap(inds[0], inds[1]);
-    }
-
-    void apply(const tc::Cosets &table, int gen) {
-        for (auto &ind : inds) {
-            ind = table.get(ind, gen);
-        }
-        flip();
-    }
-};
-
-/**
- * Produce a list of all generators for the group context. The range [0..group.ngens).
- */
-std::vector<int> generators(const tc::Group &context) {
-    std::vector<int> g_gens(context.ngens);
-    std::iota(g_gens.begin(), g_gens.end(), 0);
-    return g_gens;
-}
-
-/**
- * Determine which of g_gens are the correct names for sg_gens within the current context
- */
-std::vector<int> recontext_gens(
-    const tc::Group &context,
-    std::vector<int> g_gens,
-    std::vector<int> sg_gens) {
-
-    std::sort(g_gens.begin(), g_gens.end());
-
-    int inv_gen_map[context.ngens];
-    for (size_t i = 0; i < g_gens.size(); i++) {
-        inv_gen_map[g_gens[i]] = i;
-    }
-
-    std::vector<int> s_sg_gens;
-    s_sg_gens.reserve(sg_gens.size());
-    for (const auto gen : sg_gens) {
-        s_sg_gens.push_back(inv_gen_map[gen]);
-    }
-    std::sort(s_sg_gens.begin(), s_sg_gens.end());
-
-    return s_sg_gens;
-}
-
-/**
- * Determine whether the orientation of the group sg_gens is reversed from the group g_gens within group context
- */
-int get_parity(
-    const tc::Group &context,
-    const std::vector<int> &g_gens,
-    const std::vector<int> &sg_gens
-) {
-    if (g_gens.size() != sg_gens.size() + 1) return 0;
-
-    const auto proper_sg_gens = recontext_gens(context, g_gens, sg_gens);
-
-    int i = 0;
-    for (; i < sg_gens.size(); ++i) {
-        if (proper_sg_gens[i] != i) {
-            break;
-        }
-    }
-
-    return i & 1;
-}
-
-/**
- * Solve the cosets generated by sg_gens within the subgroup generated by g_gens of the group context
- */
-tc::Cosets solve(
-    const tc::Group &context,
-    const std::vector<int> &g_gens,
-    const std::vector<int> &sg_gens
-) {
-    const auto proper_sg_gens = recontext_gens(context, g_gens, sg_gens);
-    return context.subgroup(g_gens).solve(proper_sg_gens);
-}
-
-/**
- * Apply some context transformation to all primitives of this mesh.
- */
-template<unsigned N>
-std::vector<Primitive<N>> apply(std::vector<Primitive<N>> prims, const tc::Cosets &table, int gen) {
-    for (auto &prim : prims) {
-        prim.apply(table, gen);
-    }
-    return prims;
-}
-
-/**
- * Reverse the orientation of all primitives in this mesh.
- */
-template<unsigned N>
-void flip(std::vector<Primitive<N>> prims) {
-    for (auto &prim : prims) {
-        prim.flip();
-    }
-}
-
-/**
- * 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<unsigned N>
-[[nodiscard]]
-std::vector<Primitive<N>> recontext(
-    std::vector<Primitive<N>> prims,
-    const tc::Group &context,
-    const std::vector<int> &g_gens,
-    const std::vector<int> &sg_gens
-) {
-    const auto proper_sg_gens = recontext_gens(context, g_gens, sg_gens);
-    const auto table = solve(context, g_gens, {});
-    const auto path = solve(context, sg_gens, {}).path;
-
-    auto map = path.template walk<int, int>(0, proper_sg_gens, [table](int coset, int gen) {
-        return table.get(coset, gen);
-    });
-
-    std::vector<Primitive<N>> res(prims);
-    for (Primitive<N> &prim : res) {
-        for (auto &ind : prim.inds) {
-            ind = map[ind];
-        }
-    }
-
-    if (get_parity(context, g_gens, sg_gens) == 1)
-        flip(res);
+using mat1 = mat<1>;
+using mat2 = mat<2>;
+using mat3 = mat<3>;
+using mat4 = mat<4>;
+using mat5 = mat<5>;
 
+mat4 orthographic(float left, float right, float bottom, float top, float front, float back) {
+    mat4 res = mat4();
+    res <<
+        2 / (right - left), 0, 0, -(right + left) / (right - left),
+            0, 2 / (top - bottom), 0, -(top + bottom) / (top - bottom),
+            0, 0, 2 / (front - back), -(front + back) / (front - back),
+            0, 0, 0, 1;
     return res;
 }
 
-/**
- * Union several meshes of the same dimension
- */
-template<unsigned N>
-std::vector<Primitive<N>> merge(const std::vector<std::vector<Primitive<N>>> &meshes) {
-    size_t size = 0;
-    for (const auto &mesh : meshes) {
-        size += mesh.size();
-    }
-
-    std::vector<Primitive<N>> res;
-    res.reserve(size);
-    for (const auto &mesh : meshes) {
-        res.insert(res.end(), mesh.begin(), mesh.end());
-    }
+mat4 perspective(float fovy, float aspect, float zNear, float zFar) {
+    float tanHalfFovy(std::tan(fovy / 2));
 
+    mat4 res = mat4::Identity();
+    res(0, 0) = 1 / (aspect * tanHalfFovy);
+    res(1, 1) = 1 / (tanHalfFovy);
+    res(2, 2) = -(zFar + zNear) / (zFar - zNear);
+    res(3, 2) = -1;
+    res(2, 3) = -(2 + zFar * zNear) / (zFar - zNear);
     return res;
 }
 
-template<unsigned N>
-[[nodiscard]]
-std::vector<std::vector<Primitive<N>>> each_tile(
-    std::vector<Primitive<N>> prims,
-    const tc::Group &context,
-    const std::vector<int> &g_gens,
-    const std::vector<int> &sg_gens
-) {
-    std::vector<Primitive<N>> base = recontext(prims, context, g_gens, sg_gens);
-    const auto proper_sg_gens = recontext_gens(context, g_gens, sg_gens);
-
-    const auto table = solve(context, g_gens, {});
-    const auto path = solve(context, g_gens, sg_gens).path;
-
-    auto _gens = generators(context);
-
-    auto res = path.walk<std::vector<Primitive<N>>, int>(base, generators(context), [&](auto from, auto gen) {
-        return apply(from, table, gen);
-    });
-
+mat4 translation(float x, float y, float z) {
+    mat4 res = mat4();
+    res <<
+        1, 0, 0, x,
+        0, 1, 0, y,
+        0, 0, 1, z,
+        0, 0, 0, 1;
     return res;
 }
-
-template<unsigned N>
-[[nodiscard]]
-std::vector<Primitive<N>> tile(
-    std::vector<Primitive<N>> prims,
-    const tc::Group &context,
-    const std::vector<int> &g_gens,
-    const std::vector<int> &sg_gens
-) {
-    auto res = each_tile<N>(prims, 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<unsigned N>
-[[nodiscard]]
-std::vector<Primitive<N + 1>> fan(std::vector<Primitive<N>> prims, int root) {
-    std::vector<Primitive<N + 1>> res(prims.size());
-    std::transform(prims.begin(), prims.end(), res.begin(),
-        [root](const Primitive<N> &prim) {
-            return Primitive<N + 1>(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<unsigned N>
-std::vector<Primitive<N>> triangulate(
-    const tc::Group &context,
-    const std::vector<int> &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");
-
-    const auto &combos = Combos(g_gens, g_gens.size() - 1);
-
-    std::vector<std::vector<Primitive<N>>> meshes;
-
-    for (const auto &sg_gens : combos) {
-        auto base = triangulate<N - 1>(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));
-    }
-
-    return merge(meshes);
-}
-
-/**
- * Single-index primitives should not be further triangulated.
- */
-template<>
-std::vector<Primitive<1>> triangulate(
-    const tc::Group &context,
-    const std::vector<int> &g_gens
-) {
-    if (not g_gens.empty()) // todo make static assert
-        throw std::logic_error("g_gens must be empty for a trivial Mesh");
-
-    std::vector<Primitive<1>> res;
-    res.emplace_back();
-    return res;
-}
-
-template<unsigned N, class T>
-auto hull(const tc::Group &group, T all_sg_gens, const std::vector<std::vector<int>> &exclude) {
-    std::vector<std::vector<Primitive<N>>> parts;
-    auto g_gens = generators(group);
-    for (const std::vector<int> &sg_gens : all_sg_gens) {
-        bool excluded = false;
-        for (const auto &test : exclude) {
-            if (sg_gens == test) {
-                excluded = true;
-                break;
-            }
-        }
-        if (excluded) continue;
-
-        const auto &base = triangulate<N>(group, sg_gens);
-        const auto &tiles = each_tile(base, group, g_gens, sg_gens);
-        for (const auto &tile : tiles) {
-            parts.push_back(tile);
-        }
-    }
-    return parts;
-}

vis/include/mirror.hpp

@@ -6,189 +6,115 @@
 #include <vector>
 #include <algorithm>
 
-template<unsigned N>
-using vec = std::array<float, N>;
+#include <geometry.hpp>
 
-using vec1 = vec<1>;
-using vec2 = vec<2>;
-using vec3 = vec<3>;
-using vec4 = vec<4>;
-using vec5 = vec<5>;
+template<unsigned int N>
+Eigen::Matrix<float, N, N> mirror(const tc::Group &group) {
+    Eigen::Matrix<float, N, N> res;
+    res.setZero();
 
-template<class V>
-V operator*(V a, const float &b) {
-    for (auto &e : a) e *= b;
-    return a;
-}
+    for (int c = 0; c < group.rank(); ++c) {
+        for (int r = 0; r < c; ++r) {
+            auto angle = M_PI / group.get(c, r);
+            auto dot = res.col(c).dot(res.col(r));
 
-template<class V>
-V operator*(const float &b, V a) {
-    for (auto &e : a) e *= b;
-    return a;
-}
-
-template<class V>
-V operator/(V a, const float &b) {
-    for (auto &e : a) e /= b;
-    return a;
-}
-
-template<class V>
-V operator+(const V &a, V b) {
-    for (int i = 0; i < a.size(); ++i) {
-        a[i] += b[i];
-    }
-    return a;
-}
-
-template<class V>
-V operator-(V a, const V &b) {
-    for (int i = 0; i < a.size(); ++i) {
-        a[i] -= b[i];
-    }
-    return a;
-}
-
-template<class V>
-void operator-=(V &a, const V &b) {
-    for (int i = 0; i < a.size(); ++i) {
-        a[i] -= b[i];
-    }
-}
-
-template<class V>
-void operator+=(V &a, const V &b) {
-    for (int i = 0; i < a.size(); ++i) {
-        a[i] += b[i];
-    }
-}
-
-template<class V>
-float length(const V &a) {
-    float sum = 0;
-    for (const auto &e : a) sum += e * e;
-    return sqrtf(sum);
-}
-
-template<class V>
-V normalized(const V &a) {
-    return a / length(a);
-}
-
-template<class V>
-float dot(int n, const V &a, const V &b) {
-    float sum = 0;
-    for (int i = 0; i < n; ++i) {
-        sum += a[i] * b[i];
-    }
-    return sum;
-}
-
-template<class V>
-float dot(const V &a, const V &b) {
-    float sum = 0;
-    for (int i = 0; i < a.size(); ++i) {
-        sum += a[i] * b[i];
-    }
-    return sum;
-}
-
-template<unsigned N>
-std::vector<vec<N>> mirror(const tc::Group &group) {
-    std::vector<std::vector<float>> mirrors;
-
-    for (int p = 0; p < group.ngens; ++p) {
-        std::vector<float> vp;
-        for (int m = 0; m < p; ++m) {
-            auto &vq = mirrors[m];
-            vp.push_back((cos(M_PI / group.get(p, m)) - dot(m, vp, vq)) / vq[m]);
-        }
-        vp.push_back(std::sqrt(1 - dot(p, vp, vp)));
-
-        for (const auto &v : mirrors) {
-            if (dot(p, vp, vp) > 0) {
-                for (auto &e : vp) {
-                    e *= -1;
-                }
-                break;
-            }
+            res(r, c) = (dot - cos(angle)) / res(r, r);
         }
 
-        mirrors.push_back(vp);
+        res(c, c) = sqrt(1 - res.col(c).squaredNorm());
     }
 
-    std::vector<vec<N>> res;
-    for (const auto &v : mirrors) {
-        vec<N> rv{};
-
-        // ortho proj
-        for (int i = 0; i < std::min(v.size(), (size_t) N); ++i) {
-            rv[i] = v[i];
-        }
-
-        res.push_back(rv);
-    }
     return res;
 }
 
-template<unsigned N>
-vec<N> stereo(const vec<N + 1> &v) {
-    vec<N> r;
-    for (int i = 0; i < N; ++i) {
-        r[i] = v[i] / (1 - v[N]);
+struct Stereo {
+    template<class U>
+    auto operator()(U &&mat) const {
+        const auto Rows = std::remove_reference<U>::type::RowsAtCompileTime;
+        return std::forward<U>(mat).template topRows<Rows - 1>().rowwise() / (1 - mat.template bottomRows<1>());
     }
-    return r;
-}
+};
 
-template<class V>
-V project(const V &vec, const V &target) {
-    return dot(vec, target) / dot(target, target) * target;
-}
+struct Ortho {
+    template<class U>
+    auto operator()(U &&mat) const {
+        const auto Rows = std::remove_reference<U>::type::RowsAtCompileTime;
+        return std::forward<U>(mat).template topRows<Rows - 1>();
+    }
+};
 
-template<class V>
-V reflect(const V &a, const V &axis) {
-    return a - 2.f * project(a, axis);
-}
+struct Project {
+    template<class U, class V>
+    auto operator()(U &&point, V &&axis) const {
+        return point.dot(axis) / axis.dot(axis) * std::forward<V>(axis);
+    }
+};
 
-template<class V>
-V gram_schmidt_last(std::vector<V> vecs) {
-    for (int i = 0; i < vecs.size(); ++i) {
-        for (int j = 0; j < i; ++j) {
-            vecs[i] -= project(vecs[i], vecs[j]);
+struct Reflect {
+    template<class U, class V>
+    auto operator()(U &&point, V &&axis) const {
+        return std::forward<U>(point) - 2 * Project()(point, axis);
+    }
+};
+
+template<class Mat>
+Mat gram_schmidt(Mat mat) {
+    for (int i = 0; i < mat.cols(); ++i) {
+        for (int j = i + 1; j < mat.cols(); ++j) {
+            mat.col(j) -= Project()(mat.col(j), mat.col(i));
         }
     }
-
-    return normalized(vecs[vecs.size() - 1]);
+    return mat;
 }
 
-template<class V, class C>
-V barycentric(const std::vector<V> &basis, const C &coords) {
-    V res{};
+template<class Mat>
+Mat plane_intersections(Mat normals) {
+    auto last = normals.cols() - 1;
 
-    int N = std::min(basis.size(), coords.size());
-    for (int i = 0; i < N; ++i) {
-        res += basis[i] * coords[i];
+    Mat results(normals.rows(), normals.cols());
+    results.setZero();
+
+    Eigen::Matrix<int, Mat::ColsAtCompileTime, 1> indices(normals.cols());
+    std::iota(indices.begin(), indices.end(), 0);
+
+    for (int i = 0; i < normals.cols(); ++i) {
+        std::rotate(indices.begin(), indices.begin() + 1, indices.end());
+
+        Mat cur = normals * Eigen::PermutationMatrix<Mat::ColsAtCompileTime>(indices);
+        Mat res = gram_schmidt(cur);
+
+        results.col(i) = res.col(last);
     }
-    return res;
-}
 
-template<class V>
-std::vector<V> plane_intersections(std::vector<V> normals) {
-    std::vector<V> results(normals.size());
-
-    for (int i = 0; i < normals.size(); ++i) {
-        std::rotate(normals.begin(), normals.begin() + 1, normals.end());
-        results[i] = gram_schmidt_last(normals);
-    }
+    results.colwise().normalize();
 
     return results;
 }
 
-glm::mat4 utilRotate(const int u, const int v, const float theta) {
-    auto res = glm::identity<glm::mat4>();
-    res[u][u] = std::cos(theta);
-    res[u][v] = std::sin(theta);
-    res[v][u] = -std::sin(theta);
-    res[v][v] = std::cos(theta);
+template<unsigned N>
+mat<N> rot(int u, int v, float theta) {
+    mat<N> res = mat<N>::Identity();
+    res(u, u) = std::cos(theta);
+    res(u, v) = std::sin(theta);
+    res(v, u) = -std::sin(theta);
+    res(v, v) = std::cos(theta);
     return res;
-}
+}
+
+template<class U, class V>
+auto rotor(U const &u, V const &v) {
+    using Rotor = Eigen::Matrix<
+        typename U::Scalar,
+        std::min(U::RowsAtCompileTime, V::RowsAtCompileTime),
+        std::min(U::RowsAtCompileTime, V::RowsAtCompileTime),
+        U::Options,
+        std::max(U::MaxRowsAtCompileTime, V::MaxRowsAtCompileTime),
+        std::max(U::MaxRowsAtCompileTime, V::MaxRowsAtCompileTime)
+    >;
+
+    const auto &ident = Rotor::Identity(u.rows(), v.rows());
+    const auto &inner = (u + v) / (1 + u.dot(v)) * (u + v).transpose();
+    const auto &outer = v * u.transpose();
+
+    return ident - inner + 2 * outer;
+}

vis/include/solver.hpp

@@ -0,0 +1,192 @@
+#pragma once
+
+#include <tc/core.hpp>
+#include <cmath>
+#include <optional>
+#include <numeric>
+#include <iostream>
+
+#include <geometry.hpp>
+
+/**
+ * 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<int N>
+[[nodiscard]]
+Indices<N> recontext(
+    Indices<N> prims,
+    const tc::Group &G,
+    const tc::Group &H
+) {
+    const auto table = G.solve({});
+    const auto cosets = H.solve({});
+
+    tc::Path<size_t> path(cosets, H.gens());
+
+    std::vector<size_t> map(path.order());
+    path.walk(0, [&table](size_t coset, size_t gen) {
+        return table.get(coset, gen);
+    }, map.begin());
+
+    Indices<N> res(prims);
+    auto data = res.data();
+    for (int i = 0; i < prims.size(); ++i) {
+        data[i] = map[data[i]];
+    }
+
+    return res;
+}
+
+/**
+ * Union several meshes of the same dimension
+ */
+template<int N>
+Indices<N> merge(
+    const std::vector<Indices<N>> &meshes
+) {
+    size_t cols = 0;
+    for (const auto &mesh: meshes) {
+        cols += mesh.cols();
+    }
+
+    Indices<N> res(N, cols);
+
+    size_t offset = 0;
+    for (const Indices<N> &mesh: meshes) {
+        res.middleCols(offset, mesh.cols()) = mesh;
+        offset += mesh.cols();
+    }
+
+    return res;
+}
+
+template<int N>
+[[nodiscard]]
+std::vector<Indices<N>> tile(
+    Indices<N> prims,
+    const tc::Group &G,
+    const tc::Group &H
+) {
+    const auto &table = G.solve({});
+    const auto &cosets = G.solve(H.gens());
+
+    tc::Path<> path(cosets);
+
+    std::vector<Indices<N>> res(path.order());
+
+    path.walk(prims, [&](Indices<N> from, auto gen) {
+        for (int i = 0; i < from.size(); ++i) {
+            from(i) = table.get(from(i), gen);
+        }
+        return from;
+    }, res.begin());
+
+    return res;
+}
+
+/**
+ * Produce a mesh of higher dimension by fanning a single point to all primitives in this mesh.
+ */
+template<int N>
+[[nodiscard]]
+Indices<N + 1> fan(
+    Indices<N> prims,
+    int root
+) {
+    Indices<N + 1> res(N + 1, prims.cols());
+
+    res.topRows(N) = prims;
+    res.bottomRows(1) = 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<int N>
+Indices<N> cell(
+    const tc::Group &G
+) {
+    assert(G.rank() + 1 == N);  // g_gens size must be one less than N
+
+    std::vector<Indices<N - 1>> facets;
+
+    for (auto H: G.subs(N - 2)) {
+        auto sub_base = cell<N - 1>(H);
+        auto base = recontext(sub_base, G, H);
+        auto tiles = tile(base, G, H);
+
+        facets.insert(
+            facets.end(),
+            tiles.begin() + 1, // skip "near" facet
+            tiles.end()
+        );
+    }
+
+    return fan(merge(facets), 0);
+}
+
+/**
+ * Single-index primitives should not be further triangulated.
+ */
+template<>
+Indices<1> cell<1>(
+    const tc::Group &G
+) {
+    assert(G.rank() == 0); // rank must be 0 for trivial Mesh
+
+    return Indices<1>::Zero(1, 1);
+}
+
+template<int N>
+struct Hull {
+    struct Tiling {
+        Eigen::Index first;
+        Eigen::Index count;
+    };
+
+    std::vector<Tiling> tilings{};
+    std::vector<std::vector<size_t>> subgroups{};
+
+    Indices<N> inds{};
+
+    explicit Hull(tc::Group const &G) {
+        std::vector<Indices<N>> parts;
+
+        Eigen::Index first = 0;
+        for (const auto &H: G.subs(N - 1)) {
+            auto sub_base = cell<N>(H);
+            auto base = recontext(sub_base, G, H);
+            auto tiles = tile(base, G, H);
+            auto tiling = merge(tiles);
+
+            subgroups.push_back(H.gens());
+            tilings.push_back({first, tiling.cols()});
+
+            parts.push_back(tiling);
+            first += tiling.cols();
+        }
+
+        inds = merge(parts);
+    }
+};
+
+struct Points {
+    Eigen::Array<float, 4, Eigen::Dynamic> verts;
+
+    explicit Points(tc::Group const &G, Eigen::Vector<float, 5> root) {
+        auto mirrors = mirror<5>(G);
+        auto corners = plane_intersections(mirrors);
+
+        auto start = corners * root;
+
+        tc::Cosets table = G.solve();
+        tc::Path<vec5> path(table, mirrors.colwise());
+
+        Eigen::Array<float, 5, Eigen::Dynamic> higher(5, path.order());
+        path.walk(start, Reflect(), higher.matrix().colwise().begin());
+
+        verts = Stereo()(higher);
+    }
+};

vis/shaders/CMakeLists.txt

@@ -0,0 +1,5 @@
+add_embed_library(shaders
+    solid.fs.glsl
+    deferred.vs.glsl
+    slice.gm.glsl
+    )

vis/shaders/curve-ortho.gm.glsl

@@ -1,38 +0,0 @@
-#version 430
-
-#define SUBS 20
-
-layout(lines) in;
-layout(line_strip, max_vertices = SUBS) out;
-
-layout(std430, binding=1) buffer Positions {
-    vec4 verts[];
-};
-
-layout(std140, binding=1) uniform Matrices {
-    mat4 proj;
-    mat4 view;
-};
-
-out gl_PerVertex {
-    vec4 gl_Position;
-    float gl_PointSize;
-};
-
-in vec4 gpos[];
-
-out vec4 vpos;
-
-vec4 stereo(vec4 v) {
-    return vec4(v.xyz, 1);
-}
-
-void main() {
-    for (int i = 0; i < SUBS; i++) {
-        vpos = mix(gpos[0], gpos[1], i * 1.0f / (SUBS - 1));
-        vpos = normalize(vpos);
-        gl_Position = proj * stereo(vpos);
-        EmitVertex();
-    }
-    EndPrimitive();
-}

vis/shaders/curve-stereo.gm.glsl

@@ -1,39 +0,0 @@
-#version 430
-
-#define SUBS 20
-
-layout(lines) in;
-layout(line_strip, max_vertices = SUBS) out;
-
-layout(std430, binding=1) buffer Positions {
-    vec4 verts[];
-};
-
-layout(std140, binding=1) uniform Matrices {
-    mat4 proj;
-    mat4 view;
-};
-
-out gl_PerVertex {
-    vec4 gl_Position;
-    float gl_PointSize;
-};
-
-in vec4 gpos[];
-
-out vec4 vpos;
-
-vec4 stereo(vec4 v) {
-    return vec4(v.xyz / (1 - v.w), 1);
-}
-
-void main() {
-    for (int i = 0; i < SUBS; i++) {
-        vpos = mix(gpos[0], gpos[1], i * 1.0f / (SUBS - 1));
-        vpos = normalize(vpos);
-
-        gl_Position = proj * stereo(vpos);
-        EmitVertex();
-    }
-    EndPrimitive();
-}

vis/shaders/deferred.vs.glsl

@@ -9,11 +9,16 @@ layout(std140, binding=1) uniform Matrices {
     mat4 view;
 };
 
+layout(std140, binding=2) uniform ModelMatrices {
+    mat4 linear;
+    vec4 translation;
+};
+
 layout(location=0) in ivec4 inds;
-layout(location=1) in vec4 col;
+layout(location=1) in vec3 col;
 
 out ivec4 vInds;
-out vec4 vCol;
+out vec3 vCol;
 
 out gl_PerVertex {
     vec4 gl_Position;
@@ -24,7 +29,7 @@ void main() {
     vInds = inds;
     vCol = col;
 
-    vec4 pos = view * verts[vInds.x];
+    vec4 pos = linear * verts[vInds.x] + translation;
     gl_Position = proj * vec4(pos.xyz, 1);
     gl_PointSize = 5;
 }

vis/shaders/diffuse.fs.glsl

@@ -1,14 +0,0 @@
-#version 430
-
-layout(location=1) in vec4 col;
-layout(location=2) in vec3 normal;
-
-out vec4 color;
-
-void main() {
-    float bright = dot(normal, normalize(vec3(-0.6, 1, 2)));
-    bright = .6 + .3 * bright;
-
-    color = col;
-    color.xyz *= bright;
-}

vis/shaders/direct-ortho.vs.glsl

@@ -1,24 +0,0 @@
-#version 430
-
-layout(std430, binding=1) buffer Positions {
-    vec4 verts[];
-};
-
-layout(std140, binding=1) uniform Matrices {
-    mat4 proj;
-    mat4 view;
-};
-
-out gl_PerVertex {
-    vec4 gl_Position;
-    float gl_PointSize;
-};
-
-out vec4 vpos;
-
-void main() {
-    vec4 pos = verts[gl_VertexID];
-    vpos = view * pos;
-    gl_Position = proj * vec4(vpos.xyz / (1), 1);
-    gl_PointSize = 5;
-}

vis/shaders/direct-stereo.vs.glsl

@@ -1,24 +0,0 @@
-#version 430
-
-layout(std430, binding=1) buffer Positions {
-    vec4 verts[];
-};
-
-layout(std140, binding=1) uniform Matrices {
-    mat4 proj;
-    mat4 view;
-};
-
-out gl_PerVertex {
-    vec4 gl_Position;
-    float gl_PointSize;
-};
-
-out vec4 vpos;
-
-void main() {
-    vec4 pos = verts[gl_VertexID];
-    vpos = view * pos;
-    gl_Position = proj * vec4(vpos.xyz / (1 - vpos.w), 1);
-    gl_PointSize = 5;
-}

vis/shaders/slice.gm.glsl

@@ -12,11 +12,16 @@ layout(std140, binding=1) uniform Matrices {
     mat4 view;
 };
 
+layout(std140, binding=2) uniform ModelMatrices {
+    mat4 linear;
+    vec4 translation;
+};
+
 in ivec4 vInds[];
-in vec4 vCol[];
+in vec3 vCol[];
 
 layout(location=0) out vec4 pos;
-layout(location=1) out vec4 col;
+layout(location=1) out vec3 col;
 layout(location=2) out vec3 normal;
 
 out gl_PerVertex {
@@ -36,7 +41,8 @@ void emit(vec4 v) {
 
 void main() {
     vec4 pos4[4];
-    for (int i = 0; i < 4; ++i) pos4[i] = view * verts[vInds[0][i]];
+    for (int i = 0; i < 4; ++i)
+        pos4[i] = view * (linear * verts[vInds[0][i]] + translation);
 
     int lo[4], L = 0;
     int hi[4], H = 0;

vis/shaders/solid.fs.glsl

@@ -1,8 +1,7 @@
 #version 430
 
-layout(location=2) uniform vec3 col;
-
 layout(location=0) in vec3 pos;
+layout(location=1) in vec3 col;
 layout(location=2) in vec3 normal;
 
 out vec4 color;

vis/src/components.hpp

@@ -0,0 +1,81 @@
+#pragma once
+
+#include <cgl/buffer.hpp>
+#include <cgl/shaderprogram.hpp>
+#include <cgl/vertexarray.hpp>
+#include <cgl/pipeline.hpp>
+
+#include <tc/groups.hpp>
+
+#include <Eigen/Eigen>
+
+#include <entt/entt.hpp>
+
+#include "mirror.hpp"
+#include "geometry.hpp"
+#include "solver.hpp"
+
+#include <shaders.hpp>
+
+namespace vis {
+    using Color = Eigen::Vector<float, 3>;
+
+    template<int R_, int D_, int G_>
+    struct Structure {
+        static constexpr auto Rank = R_;
+        static constexpr auto Dim = D_;
+        static constexpr auto Grade = G_;
+
+        using AffineDf = Eigen::Transform<float, Dim, Eigen::Affine>;
+
+        using VectorRf = Eigen::Vector<float, Rank>;
+
+        // todo cache and recompute cells/points on frame (only if changed) in a system.
+
+        tc::Group group;
+        VectorRf root;
+
+        AffineDf transform = AffineDf::Identity();
+
+        explicit Structure(tc::Group const &group, VectorRf root) :
+            group(group), root(root), transform(AffineDf::Identity()) {
+        }
+    };
+
+    template<typename Str_>
+    struct Part {
+        using Str = Str_;
+
+        entt::entity parent;
+
+        GLuint first;
+        GLuint count;
+        Color color = Color::Ones();
+        bool enabled = true;
+    };
+
+    struct Command {
+        unsigned int count, instanceCount, first, baseInstance;
+    };
+
+    template<typename Str_>
+    struct VBOs {
+        using Str = Str_;
+
+        using VectorDf = Str::VectorDf;
+        using Color = Str::Color;
+        using ArrayGui = Str::ArrayGui;
+
+        struct Uniform {
+            Eigen::Matrix4f linear;
+            Eigen::Vector4f translation;
+        };
+
+        cgl::Buffer<VectorDf> vertices;
+        cgl::Buffer<Color> colors;
+        cgl::Buffer<ArrayGui> indices;
+        cgl::Buffer<Uniform> uniform;
+        cgl::Buffer<Command> commands;
+    };
+
+}

vis/src/main.cpp

@@ -1,23 +1,29 @@
 #include <glad/glad.h>
 #include <GLFW/glfw3.h>
+
 #include <cmath>
 #include <iostream>
 
-#include <glm/gtc/type_ptr.hpp>
+#include <imgui.h>
+#include <backends/imgui_impl_glfw.h>
+#include <backends/imgui_impl_opengl3.h>
 
 #include <tc/groups.hpp>
 
 #include "util.hpp"
 #include "mirror.hpp"
-#include "geometry.hpp"
 
-#include <cgl/vertexarray.hpp>
-#include <cgl/shaderprogram.hpp>
-#include <cgl/pipeline.hpp>
-#include <random>
+#include "components.hpp"
+#include "fmt/core.h"
+#include "fmt/ranges.h"
+#include "systems.hpp"
 
-#include <chrono>
-#include <yaml-cpp/yaml.h>
+#include <shaders.hpp>
+
+#ifndef NDEBUG
+#include <cgl/debug.hpp>
+#include <utility>
+#endif
 
 #ifdef _WIN32
 extern "C" {
@@ -25,394 +31,167 @@ __attribute__((unused)) __declspec(dllexport) int NvOptimusEnablement = 0x000000
 }
 #endif
 
-struct Matrices {
-    glm::mat4 proj;
-    glm::mat4 view;
+#ifndef M_PI_2f32
+#define M_PI_2f32 3.14159265358979323846f
+#endif
 
-    Matrices(const glm::mat4 &proj, const glm::mat4 &view)
-        : proj(proj), view(view) {
+struct Matrices {
+    Eigen::Matrix4f proj;
+    Eigen::Matrix4f view;
+
+    Matrices(Eigen::Matrix4f proj, Eigen::Matrix4f view)
+        : proj(std::move(proj)), view(std::move(view)) {
     }
 };
 
 struct State {
-    float time;
-    float time_delta;
-
-    float st;
-
     int dimension;
+
+    float time = 0;
+
+    Eigen::Matrix4f view = Eigen::Matrix4f::Identity();
+
+    entt::registry registry;
 };
 
-Matrices build(GLFWwindow *window, State &state) {
+Matrices build(GLFWwindow* window, State &state, ImGuiContext* ctx) {
     int width, height;
     glfwGetFramebufferSize(window, &width, &height);
 
     auto aspect = (float) width / (float) height;
     auto pheight = 1.4f;
     auto pwidth = aspect * pheight;
-    glm::mat4 proj = glm::ortho(-pwidth, pwidth, -pheight, pheight, -10.0f, 10.0f);
+    Eigen::Matrix4f proj = orthographic(-pwidth, pwidth, -pheight, pheight, -10.0f, 10.0f);
 
-    if (!glfwGetKey(window, GLFW_KEY_LEFT_SHIFT)) {
-        state.st += state.time_delta / 8;
+    auto &io = ImGui::GetIO();
+
+    if (io.MouseDown[0] && !io.WantCaptureMouse) {
+        Eigen::Vector4f src{0, 0, 1, 0};
+        Eigen::Vector4f dst{io.MouseDelta.x, -io.MouseDelta.y, 300, 0};
+        dst.normalize();
+
+        auto rotate = rotor(src, dst);
+        state.view = rotate * state.view;
     }
 
-    auto view = glm::identity<glm::mat4>();
-    if (state.dimension < 4) {
-        view *= utilRotate(2, 3, M_PI_2f32 + 0.01f);
-    }
-
-    if (state.dimension > 1) {
-        view *= utilRotate(0, 1, state.st * .40f);
-    }
-    if (state.dimension > 2) {
-        view *= utilRotate(0, 2, state.st * .20f);
-        view *= utilRotate(1, 2, state.st * .50f);
-    }
-    if (state.dimension > 3) {
-        view *= utilRotate(0, 3, state.st * 1.30f);
-        view *= utilRotate(1, 3, state.st * .25f);
-        view *= utilRotate(2, 3, state.st * 1.42f);
-    }
-
-    return Matrices(proj, view);
+    return Matrices(proj, state.view);
 }
 
-template<class C>
-std::vector<vec4> points(const tc::Group &group, const C &coords) {
-    auto cosets = group.solve();
-    auto mirrors = mirror<5>(group);
+void show_overlay(State &state) {
+    static std::string gl_vendor = (const char*) glGetString(GL_VENDOR);
+    static std::string gl_renderer = (const char*) glGetString(GL_RENDERER);
+    static std::string gl_version = (const char*) glGetString(GL_VERSION);
+    static std::string glsl_version = (const char*) glGetString(GL_SHADING_LANGUAGE_VERSION);
 
-    auto corners = plane_intersections(mirrors);
+    ImGuiWindowFlags window_flags =
+        ImGuiWindowFlags_AlwaysAutoResize |
+        ImGuiWindowFlags_NoSavedSettings |
+        ImGuiWindowFlags_NoFocusOnAppearing |
+        ImGuiWindowFlags_NoNav |
+        ImGuiWindowFlags_NoBringToFrontOnFocus |
+        ImGuiWindowFlags_NoMove;
 
-    auto start = barycentric(corners, coords);
+    ImGuiStyle &style = ImGui::GetStyle();
+    const auto PAD = style.DisplaySafeAreaPadding;
+    auto window_pos = PAD;
 
-    const auto &higher = cosets.path.walk<vec5, vec5>(start, mirrors, reflect<vec5>);
-    std::vector<vec4> lower(higher.size());
-    std::transform(higher.begin(), higher.end(), lower.begin(), stereo<4>);
-    return lower;
+    ImGui::SetNextWindowPos(window_pos, ImGuiCond_Always);
+    ImGui::SetNextWindowBgAlpha(0.35f * style.Alpha);
+    ImGui::SetNextWindowCollapsed(true, ImGuiCond_Appearing);
+    ImGui::Begin("Graphics Information", nullptr, window_flags);
+    ImGui::Text("GL Vendor    | %s", gl_vendor.c_str());
+    ImGui::Text("GL Renderer  | %s", gl_renderer.c_str());
+    ImGui::Text("GL Version   | %s", gl_version.c_str());
+    ImGui::Text("GLSL Version | %s", glsl_version.c_str());
+
+    auto v2 = ImGui::GetWindowSize();
+    window_pos.y += v2.y + PAD.y;
+    ImGui::End();
 }
 
-template<unsigned N>
-struct Prop {
-    cgl::VertexArray vao;
-    cgl::Buffer<vec4> vbo;
-    cgl::Buffer<Primitive<N>> ibo;
+void show_options(entt::registry &registry) {
+    using Slice = vis::Structure<5, 4, 4>;
+    auto view = registry.view<Slice>();
 
-    Prop() : vao(), vbo(), ibo() {}
-};
+    for (auto [entity, structure]: view.each()) {
+        ImGui::Begin("Structure View Options");
 
-template<unsigned N>
-struct Renderer {
-    std::vector<Prop<N>> props;
+        for (int idx = 0; idx < structure.parts.size(); ++idx) {
+            const auto &part_entity = structure.parts[idx];
+            auto &part = registry.get<Slice::Part>(part_entity);
 
-    virtual void bound(const std::function<void()> &action) const = 0;
+            std::string label = fmt::format("{}", idx);
 
-    virtual void _draw(const Prop<N> &) const = 0;
+            ImGui::Checkbox(label.c_str(), (bool*) (&(part.enabled)));
+            ImGui::ColorEdit3(label.c_str(), part.color.data(), ImGuiColorEditFlags_NoLabel);
+        }
 
-    void render() const {
-        bound([&]() {
-            for (const auto &prop : props) {
-                _draw(prop);
-            }
-        });
+        ImGui::End();
     }
-};
+}
 
-template<unsigned N>
-struct SliceProp : public Prop<N> {
-    vec3 color;
+void set_style() {
+    ImGui::StyleColorsDark();
 
-    SliceProp(vec3 color) : Prop<N>(), color(color) {}
+    ImGuiStyle &style = ImGui::GetStyle();
+    style.WindowRounding = 4;
+    style.FrameRounding = 2;
+    style.DisplaySafeAreaPadding.x = 10;
+    style.DisplaySafeAreaPadding.y = 10;
+}
 
-    SliceProp(SliceProp &) = delete;
-
-    SliceProp(SliceProp &&) noexcept = default;
-
-    template<class T, class C>
-    static SliceProp<N> build(
-        const tc::Group &g,
-        const C &coords,
-        vec3 color,
-        T all_sg_gens,
-        const std::vector<std::vector<int>> &exclude
-    ) {
-        SliceProp<N> res(color);
-
-        res.vbo.put(points(g, coords));
-        res.ibo.put(merge<N>(hull<N>(g, all_sg_gens, exclude)));
-        res.vao.ipointer(0, res.ibo, 4, GL_UNSIGNED_INT);
-
-        return res;
-    }
-};
-
-template<unsigned N>
-struct SliceRenderer : public Renderer<N> {
-    cgl::pgm::vert defer = cgl::pgm::vert::file(
-        "shaders/slice/deferred.vs.glsl");
-    cgl::pgm::geom slice = cgl::pgm::geom::file(
-        "shaders/slice/slice.gm.glsl");
-    cgl::pgm::frag solid = cgl::pgm::frag::file(
-        "shaders/solid.fs.glsl");
-
-    cgl::pipeline pipe;
-
-    SliceRenderer() {
-        pipe.stage(defer);
-        pipe.stage(slice);
-        pipe.stage(solid);
+int run(GLFWwindow* window, ImGuiContext* ctx) {
+    if (glfwRawMouseMotionSupported()) {
+        glfwSetInputMode(window, GLFW_RAW_MOUSE_MOTION, GLFW_TRUE);
     }
 
-    SliceRenderer(SliceRenderer &) = delete;
-
-    SliceRenderer(SliceRenderer &&) noexcept = default;
-
-    void bound(const std::function<void()> &action) const override {
-        pipe.bound(action);
-    }
-
-    void _draw(const Prop<N> &prop) const override {
-        glBindBufferBase(GL_SHADER_STORAGE_BUFFER, 1, prop.vbo);
-//        glProgramUniform3fv(solid, 2, 1, &prop.color.front());
-        glProgramUniform3f(solid, 2, 1.f, 1.f, 1.f);
-        prop.vao.bound([&]() {
-            glDrawArrays(GL_POINTS, 0, prop.ibo.count() * N);
-        });
-    }
-};
-
-
-template<unsigned N>
-struct DirectRenderer : public Renderer<N> {
-    cgl::pipeline pipe;
-
-    DirectRenderer() = default;
-
-    DirectRenderer(DirectRenderer &) = delete;
-
-    DirectRenderer(DirectRenderer &&) noexcept = default;
-
-    void bound(const std::function<void()> &action) const override {
-        pipe.bound(action);
-    }
-
-    void _draw(const Prop<N> &prop) const override {
-        glBindBufferBase(GL_SHADER_STORAGE_BUFFER, 1, prop.vbo);
-//        glProgramUniform3fv(sh.solid, 2, 1, &wire.color.front());
-        prop.vao.bound([&]() {
-            prop.ibo.bound(GL_ELEMENT_ARRAY_BUFFER, [&]() {
-                glDrawElements(GL_LINES, prop.ibo.count() * N, GL_UNSIGNED_INT, nullptr);
-            });
-        });
-    }
-};
-
-//struct WireframeRenderer : public DirectRenderer<2> {
-//    WireframeRenderer() : DirectRenderer<2>() {
-//        cgl::pgm::vert direct_stereo = cgl::pgm::vert::file(
-//            "shaders/direct-ortho.vs.glsl");
-//        cgl::pgm::frag solid = cgl::pgm::frag::file(
-//            "shaders/solid.fs.glsl");
-//
-//        glProgramUniform3f(solid, 2, .3f, .3f, .3f);
-//
-//        this->pipe.stage(direct_stereo);
-//        this->pipe.stage(solid);
-//    }
-//};
-//
-//struct WireframeStereoRenderer : public WireframeRenderer {
-//    WireframeStereoRenderer() : WireframeRenderer() {
-//        cgl::pgm::vert direct_stereo = cgl::pgm::vert::file(
-//            "shaders/direct-stereo.vs.glsl");
-//        cgl::pgm::frag solid = cgl::pgm::frag::file(
-//            "shaders/solid.fs.glsl");
-//
-//        glProgramUniform3f(solid, 2, .3f, .3f, .4f);
-//
-//        this->pipe.stage(direct_stereo);
-//        this->pipe.stage(solid);
-//    }
-//};
-//
-//struct WireframeStereoCurveRenderer : public WireframeStereoRenderer {
-//    WireframeStereoCurveRenderer() : WireframeStereoRenderer() {
-//        cgl::pgm::vert direct_stereo = cgl::pgm::vert::file(
-//            "shaders/direct-stereo.vs.glsl");
-//        cgl::pgm::geom curve = cgl::pgm::geom::file(
-//            "shaders/curve-stereo.gm.glsl"
-//        );
-//        cgl::pgm::frag solid = cgl::pgm::frag::file(
-//            "shaders/solid.fs.glsl");
-//
-//        glProgramUniform3f(solid, 2, .4f, .3f, .3f);
-//
-//        this->pipe.stage(direct_stereo);
-//        this->pipe.stage(curve);
-//        this->pipe.stage(solid);
-//    }
-//};
-
-struct WireframeProp : public Prop<2> {
-    vec3 color;
-
-    WireframeProp(vec3 color) : Prop<2>(), color(color) {}
-
-    WireframeProp(WireframeProp &) = delete;
-
-    WireframeProp(WireframeProp &&) noexcept = default;
-
-    template<class T, class C>
-    static WireframeProp build(const tc::Group &g,
-        const C &coords,
-        bool curve,
-        bool ortho,
-        vec3 color,
-        T all_sg_gens,
-        const std::vector<std::vector<int>> &exclude
-    ) {
-        WireframeProp res(color);
-
-        res.vbo.put(points(g, coords));
-        res.ibo.put(merge<2>(hull<2>(g, all_sg_gens, exclude)));
-
-        return res;
-    }
-};
-
-void run(const std::string &config_file, GLFWwindow *window) {
     glEnable(GL_DEPTH_TEST);
 
     glEnable(GL_BLEND);
     glBlendFunc(GL_SRC_ALPHA, GL_ONE_MINUS_SRC_ALPHA);
 
-    SliceRenderer<4> sRen{};
-
-    cgl::pgm::vert o = cgl::pgm::vert::file("shaders/direct-ortho.vs.glsl");
-    cgl::pgm::vert s = cgl::pgm::vert::file("shaders/direct-stereo.vs.glsl");
-    cgl::pgm::geom co = cgl::pgm::geom::file("shaders/curve-ortho.gm.glsl");
-    cgl::pgm::geom cs = cgl::pgm::geom::file("shaders/curve-stereo.gm.glsl");
-    cgl::pgm::frag solid = cgl::pgm::frag::file("shaders/solid.fs.glsl");
-    glProgramUniform3f(solid, 2, 1.f, 0.f, 0.f);
-
-    DirectRenderer<2> woRen{};
-    woRen.pipe.stage(o);
-    woRen.pipe.stage(solid);
-
-    DirectRenderer<2> wocRen{};
-    wocRen.pipe.stage(o);
-    wocRen.pipe.stage(co);
-    wocRen.pipe.stage(solid);
-
-    DirectRenderer<2> wsRen{};
-    wsRen.pipe.stage(s);
-    wsRen.pipe.stage(solid);
-
-    DirectRenderer<2> wscRen{};
-    wscRen.pipe.stage(s);
-    wscRen.pipe.stage(cs);
-    wscRen.pipe.stage(solid);
-
-    auto scene = YAML::LoadFile(config_file);
+    using Slice = vis::Structure<5, 4, 4>;
+    vis::SliceRenderer<Slice> renderer;
 
     State state{};
     glfwSetWindowUserPointer(window, &state);
 
-    state.dimension = scene["dimension"].as<int>();
+    auto &registry = state.registry;
+    state.dimension = 4;
 
-    for (const auto &group_info : scene["groups"]) {
-        auto symbol = group_info["symbol"].as<std::vector<int>>();
-        auto group = tc::schlafli(symbol);
-        auto gens = generators(group);
+    auto entity = registry.create();
+    registry.emplace<Slice>(
+        entity,
+        tc::schlafli({5, 3, 3, 2}),
+        vec5{0.80, 0.09, 0.09, 0.09, 0.09}
+    );
+    registry.emplace<vis::VBOs<Slice>>(entity);
 
-        if (group_info["slices"].IsDefined()) {
-            for (const auto &slice_info : group_info["slices"]) {
-                auto root = slice_info["root"].as<vec5>();
-                auto color = slice_info["color"].as<vec3>();
-                auto exclude = std::vector<std::vector<int>>();
-
-                if (slice_info["exclude"].IsDefined()) {
-                    exclude = slice_info["exclude"].as<std::vector<std::vector<int>>>();
-                }
-
-                if (slice_info["subgroups"].IsDefined()) {
-                    auto subgroups = slice_info["subgroups"].as<std::vector<std::vector<int>>>();
-                    sRen.props.push_back(SliceProp<4>::build(
-                        group, root, color, subgroups, exclude
-                    ));
-                } else {
-                    auto combos = Combos<int>(gens, 3);
-                    sRen.props.push_back(SliceProp<4>::build(
-                        group, root, color, combos, exclude
-                    ));
-                }
-            }
-        }
-
-        if (group_info["wires"].IsDefined()) {
-            for (const auto &wire_info : group_info["wires"]) {
-                auto root = wire_info["root"].as<vec5>();
-                auto color = wire_info["color"].as<vec3>();
-                auto exclude = std::vector<std::vector<int>>();
-                auto curve = wire_info["curve"].IsDefined() && wire_info["curve"].as<bool>();
-                auto ortho = wire_info["ortho"].IsDefined() && wire_info["ortho"].as<bool>();
-
-                if (wire_info["exclude"].IsDefined()) {
-                    exclude = wire_info["exclude"].as<std::vector<std::vector<int>>>();
-                }
-
-                if (wire_info["subgroups"].IsDefined()) {
-                    auto subgroups = wire_info["subgroups"].as<std::vector<std::vector<int>>>();
-
-                    if (ortho && curve) {
-                        wocRen.props.push_back(WireframeProp::build(
-                            group, root, curve, ortho, color, subgroups, exclude
-                        ));
-                    } else if (ortho) {
-                        woRen.props.push_back(WireframeProp::build(
-                            group, root, curve, ortho, color, subgroups, exclude
-                        ));
-                    } else if (curve) {
-                        wscRen.props.push_back(WireframeProp::build(
-                            group, root, curve, ortho, color, subgroups, exclude
-                        ));
-                    } else {
-                        wsRen.props.push_back(WireframeProp::build(
-                            group, root, curve, ortho, color, subgroups, exclude
-                        ));
-                    }
-                } else {
-                    auto combos = Combos<int>(gens, 1);
-
-                    if (ortho && curve) {
-                        wocRen.props.push_back(WireframeProp::build(
-                            group, root, curve, ortho, color, combos, exclude
-                        ));
-                    } else if (ortho) {
-                        woRen.props.push_back(WireframeProp::build(
-                            group, root, curve, ortho, color, combos, exclude
-                        ));
-                    } else if (curve) {
-                        wscRen.props.push_back(WireframeProp::build(
-                            group, root, curve, ortho, color, combos, exclude
-                        ));
-                    } else {
-                        wsRen.props.push_back(WireframeProp::build(
-                            group, root, curve, ortho, color, combos, exclude
-                        ));
-                    }
-                }
-            }
-        }
-    }
+    vis::upload_structure<Slice>(registry);
+    registry.get<Slice::Part>(
+        registry.get<Slice>(entity).parts[0]
+    ).enabled = false;
+//    registry.get<Slice>(entity).enabled[0] = false;  // disable {0,1,2} cells
 
     auto ubo = cgl::Buffer<Matrices>();
     glBindBufferBase(GL_UNIFORM_BUFFER, 1, ubo);
 
     while (!glfwWindowShouldClose(window)) {
-        auto time = (float) glfwGetTime();
-        state.time_delta = state.time - time;
-        state.time = time;
+        glfwPollEvents();
+
+        ImGui_ImplOpenGL3_NewFrame();
+        ImGui_ImplGlfw_NewFrame();
+
+        auto io = ImGui::GetIO();
+        if (io.KeysDown[GLFW_KEY_ESCAPE]) {
+            glfwSetWindowShouldClose(window, true);
+            continue;
+        }
+
+        ImGui::NewFrame();
+        show_overlay(state);
+        show_options(registry);
+        ImGui::Render();
 
         int width, height;
         glfwGetFramebufferSize(window, &width, &height);
@@ -420,26 +199,52 @@ void run(const std::string &config_file, GLFWwindow *window) {
 
         glClear(GL_COLOR_BUFFER_BIT | GL_DEPTH_BUFFER_BIT);
 
-        Matrices mats = build(window, state);
-        ubo.put(mats);
+        ubo.put(build(window, state, ctx), GL_STREAM_DRAW);
 
-        glLineWidth(1.5);
+        {
+            auto &tform = registry.get<Slice>(entity).transform;
 
-        woRen.render();
-        wsRen.render();
-        wocRen.render();
-        wscRen.render();
+            if (!io.KeysDown[GLFW_KEY_SPACE]) {
+                float speed = 1.0 / 8.0;
+                if (io.KeysDown[GLFW_KEY_LEFT_SHIFT] | io.KeysDown[GLFW_KEY_RIGHT_SHIFT]) {
+                    speed /= 4;
+                }
+                state.time += io.DeltaTime * speed;
+            }
 
-        sRen.render();
+            tform.linear().setIdentity();
+
+            if (state.dimension > 1) {
+                tform.linear() *= rot<4>(0, 1, state.time * .40f);
+            }
+            if (state.dimension > 2) {
+                tform.linear() *= rot<4>(0, 2, state.time * .20f);
+                tform.linear() *= rot<4>(1, 2, state.time * .50f);
+            }
+            if (state.dimension > 3) {
+                tform.linear() *= rot<4>(0, 3, state.time * 1.30f);
+                tform.linear() *= rot<4>(1, 3, state.time * .25f);
+                tform.linear() *= rot<4>(2, 3, state.time * 1.42f);
+            }
+
+            tform.translation().w() = std::sin(state.time * 1.4) * 1.0;
+        }
+
+        vis::upload_commands<Slice>(registry);
+        vis::upload_uniforms<Slice>(registry);
+
+        renderer(registry);
+
+        ImGui_ImplOpenGL3_RenderDrawData(ImGui::GetDrawData());
 
         glfwSwapInterval(2);
         glfwSwapBuffers(window);
-
-        glfwPollEvents();
     }
+
+    return EXIT_SUCCESS;
 }
 
-int main(int argc, char *argv[]) {
+int main(int argc, char* argv[]) {
     if (!glfwInit()) {
         std::cerr << "Failed to initialize GLFW" << std::endl;
         return EXIT_FAILURE;
@@ -449,6 +254,11 @@ int main(int argc, char *argv[]) {
     glfwWindowHint(GLFW_VERSION_MAJOR, 5);
 //    glfwWindowHint(GLFW_OPENGL_PROFILE, GLFW_OPENGL_CORE_PROFILE);
 
+    glfwWindowHint(GLFW_SAMPLES, 8);
+
+//    glfwWindowHint(GLFW_DECORATED, false);
+//    glfwWindowHint(GLFW_TRANSPARENT_FRAMEBUFFER, true);
+
     auto window = glfwCreateWindow(
         1920, 1080,
         "Coset Visualization",
@@ -463,16 +273,34 @@ int main(int argc, char *argv[]) {
     glfwMakeContextCurrent(window);
     gladLoadGLLoader((GLADloadproc) glfwGetProcAddress);
     glfwSwapInterval(1);
-    glClear(GL_COLOR_BUFFER_BIT);
-    glfwSwapBuffers(window);
 
-    std::cout << utilInfo();
+#ifndef NDEBUG
+    glEnable(GL_DEBUG_OUTPUT);
+    glDebugMessageCallback(log_gl_debug_callback, nullptr);
+#endif
 
-    std::string config_file = "presets/default.yaml";
-    if (argc > 1) config_file = std::string(argv[1]);
+    IMGUI_CHECKVERSION();
+    auto* context = ImGui::CreateContext();
+    ImGui_ImplGlfw_InitForOpenGL(window, true);
+    ImGui_ImplOpenGL3_Init("#version 130");
 
-    run(config_file, window);
+    set_style();
 
+    int exit_code = EXIT_SUCCESS;
+
+    try {
+        exit_code = run(window, context);
+    } catch (const std::exception &e) {
+        std::cerr << e.what() << std::endl;
+        exit_code = EXIT_FAILURE;
+    }
+
+    ImGui_ImplOpenGL3_Shutdown();
+    ImGui_ImplGlfw_Shutdown();
+    ImGui::DestroyContext();
+
+    glfwDestroyWindow(window);
     glfwTerminate();
-    return EXIT_SUCCESS;
+
+    return exit_code;
 }

vis/src/systems.hpp

@@ -0,0 +1,169 @@
+#pragma once
+
+#include <cgl/buffer.hpp>
+#include <cgl/shaderprogram.hpp>
+#include <cgl/vertexarray.hpp>
+#include <cgl/pipeline.hpp>
+
+#include <tc/groups.hpp>
+
+#include <Eigen/Eigen>
+
+#include <entt/entt.hpp>
+
+#include "mirror.hpp"
+#include "geometry.hpp"
+#include "solver.hpp"
+#include "components.hpp"
+
+#include <shaders.hpp>
+
+namespace vis {
+    template<typename Str>
+    void upload_structure(entt::registry &registry) {
+        {
+            auto parts = registry.view<Part<Str>>();
+            registry.destroy(parts.begin(), parts.end());
+        }
+
+        auto view = registry.view<Str, VBOs<Str>>();
+
+        for (auto [entity, structure, vbos]: view.each()) {
+            Points points(structure.group, structure.root);
+            Hull<Str::Grade> hull(structure.group);
+
+            auto &&vertices = points.verts.colwise();
+            auto &&indices = hull.inds.colwise();
+
+            vbos.vertices.put(vertices.begin(), vertices.end());
+            vbos.indices.put(indices.begin(), indices.end());
+
+            for (const auto &tiling: hull.tilings) {
+                auto part_entity = registry.create();
+                registry.emplace<Part<Str>>(
+                    part_entity,
+                    entity,
+                    tiling.first,
+                    tiling.count
+                );
+            }
+        }
+    }
+
+    template<typename Str>
+    void upload_uniforms(entt::registry &registry) {
+        auto view = registry.view<Part<Str>>();
+
+        for (auto [entity, part]: view.each()) {
+            auto &vbos = registry.get<VBOs<Str>>(part.parent);
+        }
+
+        for (auto [entity, structure, vbos]: view.each()) {
+            std::vector<typename Str::Color> colors;
+            for (const auto &part_entity: structure.parts) {
+                const auto &part = registry.get<typename Str::Part>(part_entity);
+                colors.push_back(part.color);
+            }
+            vbos.colors.put(colors.begin(), colors.end());
+
+            typename VBOs<Str>::Uniform uniform{
+                structure.transform.linear(),
+                structure.transform.translation(),
+            };
+            vbos.uniform.put(uniform, GL_STREAM_DRAW);
+        }
+    }
+
+    template<typename Str>
+    void upload_commands(entt::registry &registry) {
+        auto view = registry.view<Part<Str>>();
+        for (auto [entity, part]: view.each()) {
+            Command comm(part.count, 1, part.first, )
+        }
+
+        auto view = registry.view<Str, VBOs<Str>>();
+
+        for (auto [entity, structure, vbos]: view.each()) {
+            std::vector<typename VBOs<Str>::Command> commands;
+
+            for (unsigned int idx = 0; idx < structure.parts.size(); idx++) {
+                const auto &part_entity = structure.parts[idx];
+                const auto &part = registry.get<typename Str::Part>(part_entity);
+                if (part.enabled) {
+                    commands.push_back(
+                        {
+                            part.count,
+                            1,
+                            part.first,
+                            idx
+                        }
+                    );
+                }
+            }
+
+            vbos.commands.put(commands.begin(), commands.end(), GL_STREAM_DRAW);
+        }
+    }
+
+    template<typename Str_>
+    struct SliceRenderer {
+        using Str = Str_;
+
+        cgl::pgm::vert defer = cgl::pgm::vert(shaders::deferred_vs_glsl);
+        cgl::pgm::geom slice = cgl::pgm::geom(shaders::slice_gm_glsl);
+        cgl::pgm::frag solid = cgl::pgm::frag(shaders::solid_fs_glsl);
+
+        cgl::pipeline pipe;
+
+        cgl::VertexArray vao;
+
+        SliceRenderer() {
+            pipe.stage(defer);
+            pipe.stage(slice);
+            pipe.stage(solid);
+
+            vao.iformat(0, Str::Grade, GL_UNSIGNED_INT);  // inds
+
+            vao.format(1, 3, GL_FLOAT);  // color
+            glVertexArrayBindingDivisor(vao, 1, 1);
+        }
+
+        void operator()(entt::registry &reg) {
+            auto view = reg.view<VBOs<Str>>();
+
+            for (auto [entity, vbos]: view.each()) {
+                glBindProgramPipeline(pipe);
+
+                glBindBufferBase(
+                    GL_SHADER_STORAGE_BUFFER,
+                    1,
+                    vbos.vertices
+                );
+                glBindBufferBase(
+                    GL_UNIFORM_BUFFER,
+                    2,
+                    vbos.uniform
+                );
+
+                glBindBuffer(
+                    GL_DRAW_INDIRECT_BUFFER,
+                    vbos.commands
+                );
+
+                vao.vertexBuffer(0, vbos.indices);
+                vao.vertexBuffer(1, vbos.colors);
+
+                glBindVertexArray(vao);
+                glMultiDrawArraysIndirect(
+                    GL_POINTS,
+                    nullptr,
+                    vbos.commands.count(),
+                    0
+                );
+                glBindVertexArray(0);
+
+                glBindProgramPipeline(0);
+            }
+        }
+    };
+}