add .desktop icon to package

fix conditional threadpool include
fix windows build
2026-03-05 02:05:02 +01:00 · 2026-03-05 02:04:48 +01:00 · 2026-03-04 21:39:39 +01:00 · 2026-03-04 21:31:01 +01:00 · 2026-03-04 20:47:34 +01:00 · 2026-03-04 20:45:25 +01:00
33 changed files with 3467 additions and 1144 deletions
--- a/.gitignore
+++ b/.gitignore
@ -6,3 +6,6 @@ cmake-build-release
 /.gdb_history
 /valgrind.log
 .idea
+/perf.data
+/perf.data.old
+/clusters.puzzle
--- a/CMakeLists.txt
+++ b/CMakeLists.txt
@ -1,43 +1,51 @@
-cmake_minimum_required(VERSION 3.25)
+cmake_minimum_required(VERSION 3.28)
 project(MassSprings)

 set(CMAKE_CXX_STANDARD 26)
 set(CMAKE_EXPORT_COMPILE_COMMANDS ON)

+# Disable boost warning because our cmake/boost are recent enough
 if(POLICY CMP0167)
    cmake_policy(SET CMP0167 NEW)
 endif()

+option(DISABLE_THREADPOOL "Disable additional physics threads" OFF)
 option(DISABLE_BACKWARD "Disable backward stacktrace printer" OFF)
-option(DISABLE_TRACY "Disable the Tracy profiler client" ON)
-option(DISABLE_TESTS "Disable building and running tests" ON)
+option(DISABLE_TRACY "Disable the Tracy profiler client" OFF)
+option(DISABLE_TESTS "Disable building tests" OFF)
+option(DISABLE_BENCH "Disable building benchmarks" OFF)

-# Headers + Sources
+# Headers + Sources (excluding main.cpp)
 set(SOURCES
    src/backward.cpp
    src/graph_distances.cpp
    src/input_handler.cpp
+    src/load_save.cpp
    src/mass_spring_system.cpp
    src/octree.cpp
    src/orbit_camera.cpp
+    src/puzzle.cpp
    src/renderer.cpp
    src/state_manager.cpp
    src/threaded_physics.cpp
    src/user_interface.cpp
-    src/puzzle.cpp
 )

 # Libraries
+include(FetchContent)
 find_package(raylib REQUIRED)
-find_package(Boost REQUIRED)
-set(LIBS raylib Boost::headers)
+find_package(Boost COMPONENTS program_options REQUIRED)
+set(LIBS raylib Boost::headers Boost::program_options)
 set(FLAGS "")

 if(WIN32)
    list(APPEND LIBS opengl32 gdi32 winmm)
 endif()

-include(FetchContent)
+if(NOT DISABLE_THREADPOOL)
+    list(APPEND FLAGS THREADPOOL)
+endif()
+
 if(NOT DISABLE_BACKWARD)
    find_package(Backward REQUIRED)

@ -60,10 +68,11 @@ if(NOT DISABLE_TRACY)
    list(APPEND FLAGS TRACY)
 endif()

-# Set this after fetching tracy to hide tracy's warnings
-set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -Wall -Wextra -Wfloat-equal -Wundef -Wshadow -Wpointer-arith -Wcast-align -Wno-unused-parameter -Wunreachable-code")
+# Set this after fetching tracy to hide tracy's warnings.
+# We set -Wno-alloc-size-larger-than because it prevents BS::thread_pool from building with current gcc
+set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -Wall -Wextra -Wfloat-equal -Wundef -Wshadow -Wpointer-arith -Wcast-align -Wno-unused-parameter -Wunreachable-code -Wno-alloc-size-larger-than")
 set(CMAKE_CXX_FLAGS_DEBUG "${CMAKE_CXX_FLAGS_DEBUG} -ggdb -O0")
-set(CMAKE_CXX_FLAGS_RELEASE "${CMAKE_CXX_FLAGS_RELEASE} -ggdb -Ofast -march=native")
+set(CMAKE_CXX_FLAGS_RELEASE "${CMAKE_CXX_FLAGS_RELEASE} -ggdb -O3 -ffast-math -march=native")

 message("-- CMAKE_C_FLAGS: ${CMAKE_C_FLAGS}")
 message("-- CMAKE_C_FLAGS_DEBUG: ${CMAKE_C_FLAGS_DEBUG}")
@ -78,18 +87,21 @@ target_include_directories(masssprings PRIVATE include)
 target_link_libraries(masssprings PRIVATE ${LIBS})
 target_compile_definitions(masssprings PRIVATE ${FLAGS})

-# Testing sources
+# Testing
 if(NOT DISABLE_TESTS AND NOT WIN32)
    enable_testing()

    FetchContent_Declare(Catch2
-        GIT_REPOSITORY https://github.com/catchorg/Catch2.git
-        GIT_TAG v3.13.0
+                         GIT_REPOSITORY https://github.com/catchorg/Catch2.git
+                         GIT_TAG v3.13.0
    )
    FetchContent_MakeAvailable(Catch2)

    set(TEST_SOURCES
        test/bits.cpp
+        test/bitmap.cpp
+        test/bitmap_find_first_empty.cpp
+        # test/puzzle.cpp
    )

    add_executable(tests ${TEST_SOURCES} ${SOURCES})
@ -100,8 +112,20 @@ if(NOT DISABLE_TESTS AND NOT WIN32)
    catch_discover_tests(tests)
 endif()

+# Benchmarking
+if(NOT DISABLE_BENCH AND NOT WIN32)
+    find_package(benchmark REQUIRED)
+
+    set(BENCH_SOURCES
+        benchmark/state_space.cpp
+    )
+
+    add_executable(benchmarks ${BENCH_SOURCES} ${SOURCES})
+    target_include_directories(benchmarks PRIVATE include)
+    target_link_libraries(benchmarks benchmark raylib)
+endif()
+
 # LTO
-#if(NOT WIN32)
 include(CheckIPOSupported)
 check_ipo_supported(RESULT supported OUTPUT error)
 if(supported)
@ -110,4 +134,3 @@ if(supported)
 else()
    message(STATUS "IPO / LTO not supported")
 endif()
-#endif()
--- a/benchmark/state_space.cpp
+++ b/benchmark/state_space.cpp
@ -0,0 +1,91 @@
+#include "puzzle.hpp"
+
+#include <benchmark/benchmark.h>
+
+static std::vector<std::string> puzzles = {
+    // 0: RushHour 1
+    "S:[6x6] G:[4,2] M:[R] B:[{3x1 _ _ _ _ 1x3} {_ _ _ _ _ _} {_ _ 1x2 2X1 _ _} {_ _ _ 1x2 2x1 _} {1x2 _ 1x2 _ 2x1 _} {_ _ _ 3x1 _ _}]",
+    // 1: RushHour 2
+    "S:[6x6] G:[4,2] M:[R] B:[{1x2 3x1 _ _ 1x2 1x3} {_ 3x1 _ _ _ _} {2X1 _ 1x2 1x2 1x2 _} {2x1 _ _ _ _ _} {_ _ _ 1x2 2x1 _} {_ _ _ _ 2x1 _}]",
+    // 2: RushHour 3
+    "S:[6x6] G:[4,2] M:[R] B:[{3x1 _ _ 1x2 _ _} {1x2 2x1 _ _ _ 1x2} {_ 2X1 _ 1x2 1x2 _} {2x1 _ 1x2 _ _ 1x2} {_ _ _ 2x1 _ _} {_ 2x1 _ 2x1 _ _}]",
+    // 3: RushHour 4
+    "S:[6x6] G:[4,2] M:[R] B:[{1x3 2x1 _ _ 1x2 _} {_ 1x2 1x2 _ _ 1x3} {_ _ _ 2X1 _ _} {3x1 _ _ 1x2 _ _} {_ _ 1x2 _ 2x1 _} {2x1 _ _ 2x1 _ _}]",
+    // 4: RushHour + Walls 1
+    "S:[6x6] G:[4,2] M:[R] B:[{1x2 2x1 _ 1*1 _ _} {_ _ _ 1x2 2x1 _} {1x2 2X1 _ _ _ _} {_ _ 1x2 2x1 _ 1x3} {2x1 _ _ _ _ _} {2x1 _ 3x1 _ _ _}]",
+    // 5: RushHour + Walls 2
+    "S:[6x6] G:[4,2] M:[R] B:[{2x1 _ _ 1x2 1x2 1*1} {3x1 _ _ _ _ _} {1x2 2X1 _ 1x2 _ _} {_ _ 1x2 _ 2x1 _} {_ _ _ 2x1 _ 1x2} {_ _ 2x1 _ 1*1 _}]",
+    // 6: Dad's Puzzler
+    "S:[4x5] G:[0,3] M:[F] B:[{2X2 _ 2x1 _} {_ _ 2x1 _} {1x1 1x1 _ _} {1x2 1x2 2x1 _} {_ _ 2x1 _}]",
+    // 7: Nine Blocks
+    "S:[4x5] G:[0,3] M:[F] B:[{1x2 1x2 _ _} {_ _ 2x1 _} {1x2 1x2 2x1 _} {_ _ 2X2 _} {1x1 1x1 _ _}]",
+    // 8: Quzzle
+    "S:[4x5] G:[2,0] M:[F] B:[{2X2 _ 2x1 _} {_ _ 1x2 1x2} {_ _ _ _} {1x2 2x1 _ 1x1} {_ 2x1 _ 1x1}]",
+    // 9: Thin Klotski
+    "S:[4x5] G:[1,4] M:[F] B:[{1x2 _ 2X1 _} {_ 2x2 _ 1x1} {_ _ _ 1x1} {2x2 _ 1x1 1x1} {_ _ 1x1 1x1}]",
+    // 10: Fat Klotski
+    "S:[4x5] G:[1,3] M:[F] B:[{_ 2X2 _ 1x1} {1x1 _ _ 1x2} {1x1 2x2 _ _} {1x1 _ _ _} {1x1 1x1 2x1 _}]",
+    // 11: Klotski
+    "S:[4x5] G:[1,3] M:[F] B:[{1x2 2X2 _ 1x2} {_ _ _ _} {1x2 2x1 _ 1x2} {_ 1x1 1x1 _} {1x1 _ _ 1x1}]",
+    // 12: Century
+    "S:[4x5] G:[1,3] M:[F] B:[{1x1 2X2 _ 1x1} {1x2 _ _ 1x2} {_ 1x2 _ _} {1x1 _ _ 1x1} {2x1 _ 2x1 _}]",
+    // 13: Super Century
+    "S:[4x5] G:[1,3] M:[F] B:[{1x2 1x1 1x1 1x1} {_ 1x2 2X2 _} {1x2 _ _ _} {_ 2x1 _ 1x1} {_ 2x1 _ _}]",
+    // 14: Supercompo
+    "S:[4x5] G:[1,3] M:[F] B:[{_ 2X2 _ _} {1x1 _ _ 1x1} {1x2 2x1 _ 1x2} {_ 2x1 _ _} {1x1 2x1 _ 1x1}]",
+};
+
+static auto explore_state_space(benchmark::State& state) -> void
+{
+    const puzzle p = puzzle(puzzles[state.range(0)]);
+
+    for (auto _ : state) {
+        auto space = p.explore_state_space();
+
+        benchmark::DoNotOptimize(space);
+    }
+}
+
+static auto explore_rush_hour_puzzle_space(benchmark::State& state) -> void
+{
+    // ReSharper disable once CppTooWideScope
+    constexpr uint8_t max_blocks = 5;
+
+    constexpr uint8_t board_width = 4;
+    constexpr uint8_t board_height = 5;
+    constexpr uint8_t goal_x = board_width - 1;
+    constexpr uint8_t goal_y = 2;
+    constexpr bool restricted = true;
+
+    const boost::unordered_flat_set<puzzle::block, block_hasher2, block_equal2> permitted_blocks = {
+        puzzle::block(0, 0, 2, 1, false, false),
+        puzzle::block(0, 0, 3, 1, false, false),
+        puzzle::block(0, 0, 1, 2, false, false),
+        puzzle::block(0, 0, 1, 3, false, false)
+    };
+    const puzzle::block target_block = puzzle::block(0, 0, 2, 1, true, false);
+    constexpr std::tuple<uint8_t, uint8_t, uint8_t, uint8_t> target_block_pos_range = {
+        0,
+        goal_y,
+        board_width - 1,
+        goal_y
+    };
+
+    const puzzle p = puzzle(board_width, board_height, goal_x, goal_y, restricted, true);
+
+    for (auto _ : state) {
+        boost::unordered_flat_set<puzzle, puzzle_hasher> result = p.explore_puzzle_space(
+            permitted_blocks,
+            target_block,
+            target_block_pos_range,
+            max_blocks,
+            std::nullopt);
+
+        benchmark::DoNotOptimize(result);
+    }
+}
+
+BENCHMARK(explore_state_space)->DenseRange(0, puzzles.size() - 1)->Unit(benchmark::kMicrosecond);
+BENCHMARK(explore_rush_hour_puzzle_space)->Unit(benchmark::kSecond);
+
+BENCHMARK_MAIN();
--- a/flake.lock
+++ b/flake.lock
@ -20,11 +20,11 @@
    },
    "nixpkgs": {
      "locked": {
-        "lastModified": 1770843696,
-        "narHash": "sha256-LovWTGDwXhkfCOmbgLVA10bvsi/P8eDDpRudgk68HA8=",
+        "lastModified": 1772479524,
+        "narHash": "sha256-u7nCaNiMjqvKpE+uZz9hE7pgXXTmm5yvdtFaqzSzUQI=",
        "owner": "NixOS",
        "repo": "nixpkgs",
-        "rev": "2343bbb58f99267223bc2aac4fc9ea301a155a16",
+        "rev": "4215e62dc2cd3bc705b0a423b9719ff6be378a43",
        "type": "github"
      },
      "original": {
--- a/flake.nix
+++ b/flake.nix
@ -14,21 +14,123 @@ rec {
  # Create a shell (and possibly package) for each possible system, not only x86_64-linux
    flake-utils.lib.eachDefaultSystem (
      system: let
+        # =========================================================================================
+        # Define pkgs/stdenvs
+        # =========================================================================================
        pkgs = import nixpkgs {
          inherit system;
          config.allowUnfree = true;
          overlays = [];
        };
-        inherit (pkgs) lib stdenv;
+
+        clangPkgs = import nixpkgs {
+          inherit system;
+          config.allowUnfree = true;
+          overlays = [];
+
+          # Use this to change the compiler:
+          # - GCC: pkgs.stdenv
+          # - Clang: pkgs.clangStdenv
+          # NixOS packages are built using GCC by default. Using clang requires a full rebuild/redownload.
+          config.replaceStdenv = {pkgs}: pkgs.clangStdenv;
+        };

        windowsPkgs = import nixpkgs {
          inherit system;
+          config.allowUnfree = true;
+          overlays = [];
+
+          # Use this to cross compile to a different system
          crossSystem = {
            config = "x86_64-w64-mingw32";
          };
-          config.allowUnfree = true;
        };

+        inherit (pkgs) lib stdenv;
+
+        # =========================================================================================
+        # Define shell environment
+        # =========================================================================================
+
+        # Setup the shell when entering the "nix develop" environment (bash script).
+        shellHook = let
+          mkCmakeScript = type: let
+            typeLower = lib.toLower type;
+          in
+            pkgs.writers.writeFish "cmake-${typeLower}.fish" ''
+              cd $FLAKE_PROJECT_ROOT
+
+              echo "Removing build directory ./cmake-build-${typeLower}/"
+              rm -rf ./cmake-build-${typeLower}
+
+              echo "Creating build directory"
+              mkdir cmake-build-${typeLower}
+              cd cmake-build-${typeLower}
+
+              echo "Running cmake"
+              cmake -G "Ninja" \
+                    -DCMAKE_BUILD_TYPE="${type}" \
+                    ..
+
+              echo "Linking compile_commands.json"
+              cd ..
+              ln -sf ./cmake-build-${typeLower}/compile_commands.json ./compile_commands.json
+            '';
+
+          cmakeDebug = mkCmakeScript "Debug";
+          cmakeRelease = mkCmakeScript "Release";
+
+          mkBuildScript = type: let
+            typeLower = lib.toLower type;
+          in
+            pkgs.writers.writeFish "cmake-build.fish" ''
+              cd $FLAKE_PROJECT_ROOT/cmake-build-${typeLower}
+
+              echo "Running cmake"
+              NIX_ENFORCE_NO_NATIVE=0 cmake --build . -j$(nproc)
+            '';
+
+          buildDebug = mkBuildScript "Debug";
+          buildRelease = mkBuildScript "Release";
+
+          # Use this to specify commands that should be ran after entering fish shell
+          initProjectShell = pkgs.writers.writeFish "init-shell.fish" ''
+            echo "Entering \"${description}\" environment..."
+
+            # Determine the project root, used e.g. in cmake scripts
+            set -g -x FLAKE_PROJECT_ROOT (git rev-parse --show-toplevel)
+
+            # C/C++:
+            abbr -a cmake-debug "${cmakeDebug}"
+            abbr -a cmake-release "${cmakeRelease}"
+            abbr -a build-debug "${buildDebug}"
+            abbr -a build-release "${buildRelease}"
+            abbr -a debug-clean "${cmakeDebug} && ${buildDebug} && ./cmake-build-debug/masssprings"
+            abbr -a release-clean "${cmakeRelease} && ${buildRelease} && ./cmake-build-release/masssprings"
+            abbr -a debug "${buildDebug} && ./cmake-build-debug/masssprings"
+            abbr -a release "${buildRelease} && ./cmake-build-release/masssprings"
+
+            abbr -a run "${buildRelease} && ./cmake-build-release/masssprings"
+            abbr -a run-clusters "${buildRelease} && ./cmake-build-release/masssprings --output=clusters.puzzle --space=rh --w=6 --h=6 --gx=4 --gy=2 --blocks=4"
+            abbr -a runtests "${buildDebug} && ./cmake-build-debug/tests"
+            abbr -a runbenchs "${buildRelease} && sudo cpupower frequency-set --governor performance && ./cmake-build-release/benchmarks; sudo cpupower frequency-set --governor powersave"
+            abbr -a rungdb "${buildDebug} && gdb --tui ./cmake-build-debug/masssprings"
+            abbr -a runvalgrind "${buildDebug} && valgrind --leak-check=full --show-reachable=no --show-leak-kinds=definite,indirect,possible --track-origins=no --suppressions=valgrind.supp --log-file=valgrind.log ./cmake-build-debug/masssprings && cat valgrind.log"
+            abbr -a runperf "${buildRelease} && perf record -g ./cmake-build-release/masssprings && hotspot ./perf.data"
+            abbr -a runperf-graph "${buildRelease} && perf record -g ./cmake-build-release/benchmarks --benchmark_filter='explore_state_space' && hotspot ./perf.data"
+            abbr -a runperf-space "${buildRelease} && perf record -g ./cmake-build-release/benchmarks --benchmark_filter='explore_rush_hour_puzzle_space' && hotspot ./perf.data"
+            abbr -a runtracy "tracy -a 127.0.0.1 &; ${buildRelease} && sudo -E ./cmake-build-release/masssprings"
+
+            abbr -a runclion "clion ./CMakeLists.txt 2>/dev/null 1>&2 & disown;"
+          '';
+        in
+          builtins.concatStringsSep "\n" [
+            # Launch into pure fish shell
+            ''
+              exec "$(type -p fish)" -C "source ${initProjectShell}"
+            ''
+          ];
+
        # ===========================================================================================
        # Define custom dependencies
        # ===========================================================================================
@ -98,29 +200,20 @@ rec {
        # - Those which are needed on $PATH during the build, for example cmake and pkg-config
        # - Setup hooks, for example makeWrapper
        # - Interpreters needed by patchShebangs for build scripts (with the --build flag), which can be the case for e.g. perl
+        # NOTE: Do not add compiler here, they are provided by the stdenv
        nativeBuildInputs = with pkgs; [
          # Languages:
-          binutils
-          gcc
+          # binutils

          # C/C++:
-          gdb
-          valgrind
-          # gnumake
          cmake
          ninja
-          # cling
-          # pkg-config
-          # clang-tools
-          # compdb
-          # pprof
-          # gprof2dot
+          gdb
+          valgrind
+          kdePackages.kcachegrind
          perf
          hotspot
-          kdePackages.kcachegrind
-          gdbgui
-          massif-visualizer
-          heaptrack
+          # heaptrack
          # renderdoc
        ];

@ -134,15 +227,18 @@ rec {
          thread-pool
          boost

-          # Debugging
+          # Debugging/Testing/Profiling
          tracy-wayland
          backward-cpp
          libbfd
          catch2_3
+          gbenchmark
        ];
+
        # ===========================================================================================
        # Define buildable + installable packages
        # ===========================================================================================
+
        package = stdenv.mkDerivation rec {
          inherit buildInputs;
          pname = "masssprings";
@ -152,12 +248,17 @@ rec {
          nativeBuildInputs = with pkgs; [
            gcc
            cmake
+
+            # Fix the working directory so the auxiliary files are always available
+            makeWrapper
          ];

          cmakeFlags = [
+            "-DDISABLE_THREADPOOL=Off"
            "-DDISABLE_TRACY=On"
            "-DDISABLE_BACKWARD=On"
            "-DDISABLE_TESTS=On"
+            "-DDISABLE_BENCH=On"
          ];

          hardeningDisable = ["all"];
@ -167,11 +268,25 @@ rec {
          '';

          installPhase = ''
+            mkdir -p $out/lib
+            cp ./${pname} $out/lib/
+            cp $src/default.puzzle $out/lib/
+            cp -r $src/fonts $out/lib/fonts
+            cp -r $src/shader $out/lib/shader
+
+            # The wrapper enters the correct working dir, so fonts/shaders/presets are available
            mkdir -p $out/bin
-            cp ./${pname} $out/bin/
-            cp $src/default.puzzle $out/bin/
-            cp -r $src/fonts $out/bin/fonts
-            cp -r $src/shader $out/bin/shader
+            makeWrapper $out/lib/${pname} $out/bin/${pname} --chdir "$out/lib"
+
+            # Generate a .desktop file
+            mkdir -p $out/share/applications
+            cat <<INI > $out/share/applications/${pname}.desktop
+            [Desktop Entry]
+            Terminal=true
+            Name=PuzzleSpaces
+            Exec=$out/bin/${pname} %f
+            Type=Application
+            INI
          '';
        };

@ -190,13 +305,24 @@ rec {
            raylib
            raygui
            thread-pool
-            boost
+
+            # Disable stacktrace since that's platform dependant and won't cross compile to windows
+            # https://github.com/NixOS/nixpkgs/blob/master/pkgs/development/libraries/boost/generic.nix#L43
+            (boost.override {
+              enableShared = false;
+              extraB2Args = [
+                "--without-stacktrace"
+              ];
+            })
          ];

          cmakeFlags = [
            "-DCMAKE_SYSTEM_NAME=Windows"
+            "-DDISABLE_THREADPOOL=Off"
            "-DDISABLE_TRACY=On"
            "-DDISABLE_BACKWARD=On"
+            "-DDISABLE_TESTS=On"
+            "-DDISABLE_BENCH=On"
          ];

          installPhase = ''
@ -218,156 +344,53 @@ rec {
        # Provide environment for "nix develop"
        devShells = {
          default = pkgs.mkShell {
-            inherit nativeBuildInputs buildInputs;
+            inherit nativeBuildInputs buildInputs shellHook;
            name = description;

            # =========================================================================================
            # Define environment variables
            # =========================================================================================

-            # Custom dynamic libraries:
-            # LD_LIBRARY_PATH = builtins.concatStringsSep ":" [
-            #   # Rust Bevy GUI app:
-            #   # "${pkgs.xorg.libX11}/lib"
-            #   # "${pkgs.xorg.libXcursor}/lib"
-            #   # "${pkgs.xorg.libXrandr}/lib"
-            #   # "${pkgs.xorg.libXi}/lib"
-            #   # "${pkgs.libGL}/lib"
-            #
-            #   # JavaFX app:
-            #   # "${pkgs.libGL}/lib"
-            #   # "${pkgs.gtk3}/lib"
-            #   # "${pkgs.glib.out}/lib"
-            #   # "${pkgs.xorg.libXtst}/lib"
-            # ];
-
            # Dynamic libraries from buildinputs:
            LD_LIBRARY_PATH = nixpkgs.lib.makeLibraryPath buildInputs;
-
-            # =========================================================================================
-            # Define shell environment
-            # =========================================================================================
-
-            # Setup the shell when entering the "nix develop" environment (bash script).
-            shellHook = let
-              mkCmakeScript = type: let
-                typeLower = lib.toLower type;
-              in
-                pkgs.writers.writeFish "cmake-${typeLower}.fish" ''
-                  cd $FLAKE_PROJECT_ROOT
-
-                  # set -g -x CC ${pkgs.clang}/bin/clang
-                  # set -g -x CXX ${pkgs.clang}/bin/clang++
-
-                  echo "Removing build directory ./cmake-build-${typeLower}/"
-                  rm -rf ./cmake-build-${typeLower}
-
-                  echo "Creating build directory"
-                  mkdir cmake-build-${typeLower}
-                  cd cmake-build-${typeLower}
-
-                  echo "Running cmake"
-                  cmake -G "Ninja" \
-                        -DCMAKE_BUILD_TYPE="${type}" \
-                        ..
-
-                  echo "Linking compile_commands.json"
-                  cd ..
-                  ln -sf ./cmake-build-${typeLower}/compile_commands.json ./compile_commands.json
-                '';
-
-              cmakeDebug = mkCmakeScript "Debug";
-              cmakeRelease = mkCmakeScript "Release";
-
-              mkBuildScript = type: let
-                typeLower = lib.toLower type;
-              in
-                pkgs.writers.writeFish "cmake-build.fish" ''
-                  cd $FLAKE_PROJECT_ROOT/cmake-build-${typeLower}
-
-                  echo "Running cmake"
-                  NIX_ENFORCE_NO_NATIVE=0 cmake --build . -j$(nproc)
-                '';
-
-              buildDebug = mkBuildScript "Debug";
-              buildRelease = mkBuildScript "Release";
-
-              # Use this to specify commands that should be ran after entering fish shell
-              initProjectShell = pkgs.writers.writeFish "init-shell.fish" ''
-                echo "Entering \"${description}\" environment..."
-
-                # Determine the project root, used e.g. in cmake scripts
-                set -g -x FLAKE_PROJECT_ROOT (git rev-parse --show-toplevel)
-
-                # C/C++:
-                abbr -a cmake-debug "${cmakeDebug}"
-                abbr -a cmake-release "${cmakeRelease}"
-                abbr -a build-debug "${buildDebug}"
-                abbr -a build-release "${buildRelease}"
-                abbr -a debug "${buildDebug} && ./cmake-build-debug/masssprings"
-                abbr -a release "${buildRelease} && ./cmake-build-release/masssprings"
-                abbr -a debug-clean "${cmakeDebug} && ${buildDebug} && ./cmake-build-debug/masssprings"
-                abbr -a release-clean "${cmakeRelease} && ${buildRelease} && ./cmake-build-release/masssprings"
-
-                abbr -a rungdb "${buildDebug} && gdb --tui ./cmake-build-debug/masssprings"
-                abbr -a runperf "${buildRelease} && perf record -g ./cmake-build-release/masssprings && hotspot ./perf.data"
-                abbr -a runtracy "tracy -a 127.0.0.1 &; ${buildRelease} && sudo -E ./cmake-build-release/masssprings"
-                abbr -a runvalgrind "${buildDebug} && valgrind --leak-check=full --show-reachable=no --show-leak-kinds=definite,indirect,possible --track-origins=no --suppressions=valgrind.supp --log-file=valgrind.log ./cmake-build-debug/masssprings && cat valgrind.log"
-                abbr -a runtests "${buildDebug} && ./cmake-build-debug/tests"
-
-                abbr -a runclion "clion ./CMakeLists.txt 2>/dev/null 1>&2 & disown;"
-              '';
-            in
-              builtins.concatStringsSep "\n" [
-                # Launch into pure fish shell
-                ''
-                  exec "$(type -p fish)" -C "source ${initProjectShell} && abbr -a menu '${pkgs.bat}/bin/bat "${initProjectShell}"'"
-                ''
-              ];
          };

-          # TODO: Can't get renderdoc in FHS to work
+          # Provide environment with clang stdenv for "nix develop .#clang"
+          # TODO: Broken. Clang can't find stdlib headers or library headers (raylib, backward, ...).
+          #       Does the clangStdenv not automatically collect the include paths?
+          clang =
+            pkgs.mkShell.override {
+              stdenv = pkgs.clangStdenv;
+            } {
+              inherit shellHook;
+              name = description;

-          # FHS environment for renderdoc. Access with "nix develop .#renderdoc".
-          # https://ryantm.github.io/nixpkgs/builders/special/fhs-environments
-          # renderdoc =
-          #   (pkgs.buildFHSEnv {
-          #     name = "renderdoc-env";
-          #
-          #     targetPkgs = pkgs:
-          #       with pkgs; [
-          #         # RenderDoc
-          #         renderdoc
-          #
-          #         # Build tools
-          #         gcc
-          #         cmake
-          #
-          #         # Raylib
-          #         raylib
-          #         libGL
-          #         mesa
-          #
-          #         # X11
-          #         libx11
-          #         libxcursor
-          #         libxrandr
-          #         libxinerama
-          #         libxi
-          #         libxext
-          #         libxfixes
-          #
-          #         # Wayland
-          #         wayland
-          #         wayland-protocols
-          #         libxkbcommon
-          #       ];
-          #
-          #     runScript = "fish";
-          #
-          #     profile = ''
-          #     '';
-          #   }).env;
+              nativeBuildInputs = with pkgs; [
+                cmake
+                ninja
+              ];
+
+              buildInputs = with pkgs; [
+                # C/C++:
+                raylib
+                raygui
+                thread-pool
+                boost
+
+                # Debugging/Testing/Profiling
+                backward-cpp
+                libbfd
+                catch2_3
+                gbenchmark
+              ];
+
+              # =========================================================================================
+              # Define environment variables
+              # =========================================================================================
+
+              # Dynamic libraries from buildinputs:
+              LD_LIBRARY_PATH = nixpkgs.lib.makeLibraryPath buildInputs;
+            };
        };
      }
    );
--- a/include/config.hpp
+++ b/include/config.hpp
@ -3,12 +3,31 @@

 #include <raylib.h>

-#define THREADPOOL // Enable physics threadpool
+
+// TODO: Using the octree from the last frame completely breaks the physics :/
+// #define ASYNC_OCTREE

 // Gets set by CMake
+// #define THREADPOOL // Enable physics threadpool
 // #define BACKWARD   // Enable pretty stack traces
 // #define TRACY      // Enable tracy profiling support

+#ifdef TRACY
+#include <tracy/Tracy.hpp>
+#endif
+
+#if defined(_WIN32)
+#define NOGDI  // All GDI defines and routines
+#define NOUSER // All USER defines and routines
+#endif
+#define BS_THREAD_POOL_NATIVE_EXTENSIONS
+// ReSharper disable once CppUnusedIncludeDirective
+#include <BS_thread_pool.hpp>
+#if defined(_WIN32) // raylib uses these names as function parameters
+#undef near
+#undef far
+#endif
+
 // Window
 constexpr int INITIAL_WIDTH = 600;
 constexpr int INITIAL_HEIGHT = 600;
@ -29,7 +48,7 @@ constexpr float CAMERA_FOV = 90.0;
 constexpr float FOV_SPEED = 1.0;
 constexpr float MIN_FOV = 10.0;
 constexpr float MAX_FOV = 180.0;
-constexpr float CAMERA_DISTANCE = 20.0;
+constexpr float CAMERA_DISTANCE = 150.0;
 constexpr float ZOOM_SPEED = 2.5;
 constexpr float MIN_CAMERA_DISTANCE = 2.0;
 constexpr float MAX_CAMERA_DISTANCE = 2000.0;
--- a/include/graph_distances.hpp
+++ b/include/graph_distances.hpp
@ -21,4 +21,4 @@ public:
    [[nodiscard]] auto get_shortest_path(size_t source) const -> std::vector<size_t>;
 };

-#endif
+#endif
--- a/include/input_handler.hpp
+++ b/include/input_handler.hpp
@ -1,5 +1,5 @@
-#ifndef INPUT_HPP_
-#define INPUT_HPP_
+#ifndef INPUT_HANDLER_HPP_
+#define INPUT_HANDLER_HPP_

 #include "orbit_camera.hpp"
 #include "state_manager.hpp"
@ -78,7 +78,7 @@ public:

    // Camera
    bool camera_lock = true;
-    bool camera_mass_center_lock = false;
+    bool camera_mass_center_lock = true;
    bool camera_panning = false;
    bool camera_rotating = false;

--- a/include/load_save.hpp
+++ b/include/load_save.hpp
@ -0,0 +1,15 @@
+#ifndef LOAD_SAVE_HPP_
+#define LOAD_SAVE_HPP_
+
+#include "puzzle.hpp"
+
+#include <string>
+
+auto parse_preset_file(const std::string& preset_file) -> std::pair<std::vector<puzzle>, std::vector<std::string>>;
+auto append_preset_file(const std::string& preset_file, const std::string& preset_name, const puzzle& p) -> bool;
+auto append_preset_file_quiet(const std::string& preset_file,
+                              const std::string& preset_name,
+                              const puzzle& p,
+                              bool validate) -> bool;
+
+#endif
--- a/include/mass_spring_system.hpp
+++ b/include/mass_spring_system.hpp
@ -2,49 +2,14 @@
 #define MASS_SPRING_SYSTEM_HPP_

 #include "octree.hpp"
-#include "util.hpp"
 #include "config.hpp"

+#include <optional>
 #include <raylib.h>
-#include <raymath.h>
-
-#ifdef THREADPOOL
-#if defined(_WIN32)
-#define NOGDI  // All GDI defines and routines
-#define NOUSER // All USER defines and routines
-#endif
-#define BS_THREAD_POOL_NATIVE_EXTENSIONS
-#include <BS_thread_pool.hpp>
-#if defined(_WIN32) // raylib uses these names as function parameters
-#undef near
-#undef far
-#endif
-#endif

 class mass_spring_system
 {
 public:
-    class mass
-    {
-    public:
-        Vector3 position = Vector3Zero();
-        Vector3 previous_position = Vector3Zero(); // for verlet integration
-        Vector3 velocity = Vector3Zero();
-        Vector3 force = Vector3Zero();
-
-    public:
-        mass() = delete;
-
-        explicit mass(const Vector3 _position)
-            : position(_position), previous_position(_position) {}
-
-    public:
-        auto clear_force() -> void;
-        auto calculate_velocity(float delta_time) -> void;
-        auto calculate_position(float delta_time) -> void;
-        auto verlet_update(float delta_time) -> void;
-    };
-
    class spring
    {
    public:
@ -54,61 +19,45 @@ public:
    public:
        spring(const size_t _a, const size_t _b)
            : a(_a), b(_b) {}
-
-    public:
-        static auto calculate_spring_force(mass& _a, mass& _b) -> void;
    };

-private:
-    #ifdef THREADPOOL
-    BS::thread_pool<> threads;
-    #endif
-
 public:
+    static constexpr int SMALL_TASK_BLOCK_SIZE = 256;
+    static constexpr int LARGE_TASK_BLOCK_SIZE = 256;
+
    octree tree;

    // This is the main ownership of all the states/masses/springs.
-    std::vector<mass> masses;
+    std::vector<Vector3> positions;
+    std::vector<Vector3> previous_positions; // for verlet integration
+    std::vector<Vector3> velocities;
+    std::vector<Vector3> forces;
+
    std::vector<spring> springs;

 public:
-    mass_spring_system()
-    #ifdef THREADPOOL
-        : threads(std::thread::hardware_concurrency() - 1, set_thread_name)
-    #endif
-    {
-        infoln("Using Barnes-Hut + Octree repulsion force calculation.");
-
-        #ifdef THREADPOOL
-        infoln("Thread-pool: {} threads.", threads.get_thread_count());
-        #else
-        infoln("Thread-pool: Disabled.");
-        #endif
-    }
+    mass_spring_system() {}

    mass_spring_system(const mass_spring_system& copy) = delete;
    auto operator=(const mass_spring_system& copy) -> mass_spring_system& = delete;
    mass_spring_system(mass_spring_system& move) = delete;
    auto operator=(mass_spring_system&& move) -> mass_spring_system& = delete;

-private:
-    #ifdef THREADPOOL
-    static auto set_thread_name(size_t idx) -> void;
-    #endif
-
-    auto build_octree() -> void;
-
 public:
    auto clear() -> void;
    auto add_mass() -> void;
    auto add_spring(size_t a, size_t b) -> void;

    auto clear_forces() -> void;
-    auto calculate_spring_forces() -> void;
-    auto calculate_repulsion_forces() -> void;
-    auto verlet_update(float delta_time) -> void;
+    auto calculate_spring_force(size_t s) -> void;
+    auto calculate_spring_forces(std::optional<BS::thread_pool<>* const> thread_pool = std::nullopt) -> void;
+    auto calculate_repulsion_forces(std::optional<BS::thread_pool<>* const> thread_pool = std::nullopt) -> void;
+    auto integrate_velocity(size_t m, float dt) -> void;
+    auto integrate_position(size_t m, float dt) -> void;
+    auto verlet_update(size_t m, float dt) -> void;
+    auto update(float dt, std::optional<BS::thread_pool<>* const> thread_pool = std::nullopt) -> void;

-    auto center_masses() -> void;
+    auto center_masses(std::optional<BS::thread_pool<>* const> thread_pool = std::nullopt) -> void;
 };

 #endif
--- a/include/octree.hpp
+++ b/include/octree.hpp
@ -2,9 +2,10 @@
 #define OCTREE_HPP_

 #include <array>
+#include <vector>
+
 #include <raylib.h>
 #include <raymath.h>
-#include <vector>

 class octree
 {
@ -31,22 +32,20 @@ public:
 public:
    octree() = default;

-    octree(const octree& copy) = delete;
-    auto operator=(const octree& copy) -> octree& = delete;
-    octree(octree&& move) = delete;
-    auto operator=(octree&& move) -> octree& = delete;
+    // octree(const octree& copy) = delete;
+    // auto operator=(const octree& copy) -> octree& = delete;
+    // octree(octree&& move) = delete;
+    // auto operator=(octree&& move) -> octree& = delete;

 public:
-    auto create_empty_leaf(const Vector3& box_min, const Vector3& box_max) -> int;
-
    [[nodiscard]] auto get_octant(int node_idx, const Vector3& pos) const -> int;
-
    [[nodiscard]] auto get_child_bounds(int node_idx, int octant) const
        -> std::pair<Vector3, Vector3>;
-
+    auto create_empty_leaf(const Vector3& box_min, const Vector3& box_max) -> int;
    auto insert(int node_idx, int mass_id, const Vector3& pos, float mass, int depth) -> void;
+    static auto build_octree(octree& t, const std::vector<Vector3>& positions) -> void;

    [[nodiscard]] auto calculate_force(int node_idx, const Vector3& pos) const -> Vector3;
 };

-#endif
+#endif
--- a/include/orbit_camera.hpp
+++ b/include/orbit_camera.hpp
@ -28,4 +28,4 @@ public:
                bool mass_center_lock) -> void;
 };

-#endif
+#endif
--- a/include/puzzle.hpp
+++ b/include/puzzle.hpp
@ -1,16 +1,24 @@
 #ifndef PUZZLE_HPP_
 #define PUZZLE_HPP_

+#include "config.hpp"
 #include "util.hpp"

 #include <array>
 #include <cstddef>
-#include <format>
 #include <functional>
 #include <ranges>
 #include <string>
 #include <vector>
-#include <bits/fs_fwd.h>
+#include <boost/unordered/unordered_flat_set.hpp>
+
+// #define RUNTIME_CHECKS
+
+// Forward declare to use in puzzle member functions
+struct block_hasher;
+struct block_hasher2;
+struct block_equal2;
+struct puzzle_hasher;

 /*
 * 8x8 board
@ -125,29 +133,33 @@ public:
        }

    private:
-        [[nodiscard]] static auto create_repr(uint8_t x, uint8_t y, uint8_t w, uint8_t h, bool t = false,
-                                              bool i = false) -> uint16_t;
-
        // Repr setters
-        [[nodiscard]] auto set_x(uint8_t x) const -> block;
-        [[nodiscard]] auto set_y(uint8_t y) const -> block;
-        [[nodiscard]] auto set_width(uint8_t width) const -> block;
-        [[nodiscard]] auto set_height(uint8_t height) const -> block;
-        [[nodiscard]] auto set_target(bool target) const -> block;
-        [[nodiscard]] auto set_immovable(bool immovable) const -> block;
+        [[nodiscard]] static auto create_repr(uint8_t x,
+                                              uint8_t y,
+                                              uint8_t w,
+                                              uint8_t h,
+                                              bool t = false,
+                                              bool i = false) -> uint16_t;
+        [[nodiscard]] inline auto set_x(uint8_t x) const -> block;
+        [[nodiscard]] inline auto set_y(uint8_t y) const -> block;
+        [[nodiscard]] inline auto set_width(uint8_t width) const -> block;
+        [[nodiscard]] inline auto set_height(uint8_t height) const -> block;
+        [[nodiscard]] inline auto set_target(bool target) const -> block;
+        [[nodiscard]] inline auto set_immovable(bool immovable) const -> block;

    public:
-        [[nodiscard]] auto unpack_repr() const -> std::tuple<uint8_t, uint8_t, uint8_t, uint8_t, bool, bool>;
-
        // Repr getters
-        [[nodiscard]] auto get_x() const -> uint8_t;
-        [[nodiscard]] auto get_y() const -> uint8_t;
-        [[nodiscard]] auto get_width() const -> uint8_t;
-        [[nodiscard]] auto get_height() const -> uint8_t;
-        [[nodiscard]] auto get_target() const -> bool;
-        [[nodiscard]] auto get_immovable() const -> bool;
+        [[nodiscard]] auto unpack_repr() const -> std::tuple<uint8_t, uint8_t, uint8_t, uint8_t, bool, bool>;
+        [[nodiscard]] inline auto get_x() const -> uint8_t;
+        [[nodiscard]] inline auto get_y() const -> uint8_t;
+        [[nodiscard]] inline auto get_width() const -> uint8_t;
+        [[nodiscard]] inline auto get_height() const -> uint8_t;
+        [[nodiscard]] inline auto get_target() const -> bool;
+        [[nodiscard]] inline auto get_immovable() const -> bool;

        // Util
+        [[nodiscard]] auto hash() const -> size_t;
+        [[nodiscard]] auto position_independent_hash() const -> size_t;
        [[nodiscard]] auto valid() const -> bool;
        [[nodiscard]] auto principal_dirs() const -> uint8_t;
        [[nodiscard]] auto covers(int _x, int _y) const -> bool;
@ -199,7 +211,8 @@ private:
        // repr_cooked() = delete;
        // repr_cooked(const repr_cooked& copy) = delete;
        // repr_cooked(repr_cooked&& move) = delete;
-    } __attribute__((packed));
+    }
+        PACKED;

    /**
     * With gcc, were allowed to acces the members arbitrarily, even if they're not active (not the ones last written):
@ -247,7 +260,12 @@ public:
        }
    }

-    puzzle(const uint8_t w, const uint8_t h, const uint8_t tx, const uint8_t ty, const bool r, const bool g,
+    puzzle(const uint8_t w,
+           const uint8_t h,
+           const uint8_t tx,
+           const uint8_t ty,
+           const bool r,
+           const bool g,
           const std::array<uint16_t, MAX_BLOCKS>& b)
        : repr(create_repr(w, h, tx, ty, r, g, b))
    {
@ -259,6 +277,14 @@ public:
        }
    }

+    puzzle(const uint8_t w, const uint8_t h)
+        : repr(create_repr(w, h, 0, 0, false, false, invalid_blocks()))
+    {
+        if (w < MIN_WIDTH || w > MAX_WIDTH || h < MIN_HEIGHT || h > MAX_HEIGHT) {
+            throw std::invalid_argument("Board size out of bounds");
+        }
+    }
+
    explicit puzzle(const std::string& string_repr)
        : repr(create_repr(string_repr)) {}

@ -331,33 +357,38 @@ private:
        return blocks;
    }

-    [[nodiscard]] static auto create_meta(const std::tuple<uint8_t, uint8_t, uint8_t, uint8_t, bool, bool>& meta) -> uint16_t;
-    [[nodiscard]] static auto create_repr(uint8_t w, uint8_t h, uint8_t tx, uint8_t ty, bool r, bool g,
-                                          const std::array<uint16_t, MAX_BLOCKS>& b) -> repr_cooked;
-
-    [[nodiscard]] static auto create_repr(uint64_t byte_0, uint64_t byte_1, uint64_t byte_2,
-                                          uint64_t byte_3) -> repr_cooked;
-
-    [[nodiscard]] static auto create_repr(const std::string& string_repr) -> repr_cooked;
-
    // Repr setters
-    [[nodiscard]] auto set_restricted(bool restricted) const -> puzzle;
-    [[nodiscard]] auto set_width(uint8_t width) const -> puzzle;
-    [[nodiscard]] auto set_height(uint8_t height) const -> puzzle;
-    [[nodiscard]] auto set_goal(bool goal) const -> puzzle;
-    [[nodiscard]] auto set_goal_x(uint8_t target_x) const -> puzzle;
-    [[nodiscard]] auto set_goal_y(uint8_t target_y) const -> puzzle;
+    [[nodiscard]] static auto create_meta(
+        const std::tuple<uint8_t, uint8_t, uint8_t, uint8_t, bool, bool>& meta) -> uint16_t;
+    [[nodiscard]] static auto create_repr(uint8_t w,
+                                          uint8_t h,
+                                          uint8_t tx,
+                                          uint8_t ty,
+                                          bool r,
+                                          bool g,
+                                          const std::array<uint16_t, MAX_BLOCKS>& b) -> repr_cooked;
+    [[nodiscard]] static auto create_repr(uint64_t byte_0,
+                                          uint64_t byte_1,
+                                          uint64_t byte_2,
+                                          uint64_t byte_3) -> repr_cooked;
+    [[nodiscard]] static auto create_repr(const std::string& string_repr) -> repr_cooked;
+    [[nodiscard]] inline auto set_restricted(bool restricted) const -> puzzle;
+    [[nodiscard]] inline auto set_width(uint8_t width) const -> puzzle;
+    [[nodiscard]] inline auto set_height(uint8_t height) const -> puzzle;
+    [[nodiscard]] inline auto set_goal(bool goal) const -> puzzle;
+    [[nodiscard]] inline auto set_goal_x(uint8_t target_x) const -> puzzle;
+    [[nodiscard]] inline auto set_goal_y(uint8_t target_y) const -> puzzle;
    [[nodiscard]] auto set_blocks(std::array<uint16_t, MAX_BLOCKS> blocks) const -> puzzle;

 public:
    // Repr getters
    [[nodiscard]] auto unpack_meta() const -> std::tuple<uint8_t, uint8_t, uint8_t, uint8_t, bool, bool>;
-    [[nodiscard]] auto get_restricted() const -> bool;
-    [[nodiscard]] auto get_width() const -> uint8_t;
-    [[nodiscard]] auto get_height() const -> uint8_t;
-    [[nodiscard]] auto get_goal() const -> bool;
-    [[nodiscard]] auto get_goal_x() const -> uint8_t;
-    [[nodiscard]] auto get_goal_y() const -> uint8_t;
+    [[nodiscard]] inline auto get_restricted() const -> bool;
+    [[nodiscard]] inline auto get_width() const -> uint8_t;
+    [[nodiscard]] inline auto get_height() const -> uint8_t;
+    [[nodiscard]] inline auto get_goal() const -> bool;
+    [[nodiscard]] inline auto get_goal_x() const -> uint8_t;
+    [[nodiscard]] inline auto get_goal_y() const -> uint8_t;

    // Util
    [[nodiscard]] auto hash() const -> size_t;
@ -390,19 +421,46 @@ public:
    [[nodiscard]] auto try_move_block_at(uint8_t x, uint8_t y, direction dir) const -> std::optional<puzzle>;

    // Statespace
-    [[nodiscard]] auto try_move_block_at_fast(uint64_t bitmap, uint8_t block_idx,
-                                              direction dir) const -> std::optional<puzzle>;
-    static auto sorted_replace(std::array<uint16_t, MAX_BLOCKS> blocks, uint8_t idx,
-                               uint16_t new_val) -> std::array<uint16_t, MAX_BLOCKS>;
-    auto blocks_bitmap() const -> uint64_t;
-    static inline auto bitmap_set_bit(uint64_t bitmap, uint8_t x, uint8_t y) -> uint64_t;
-    static inline auto bitmap_get_bit(uint64_t bitmap, uint8_t x, uint8_t y) -> bool;
-    static auto bitmap_clear_block(uint64_t bitmap, block b) -> uint64_t;
-    static auto bitmap_check_collision(uint64_t bitmap, block b) -> bool;
-    static auto bitmap_check_collision(uint64_t bitmap, block b, direction dir) -> bool;
+    [[nodiscard]] INLINE inline auto try_move_block_at_fast(uint64_t bitmap,
+                                                            uint8_t block_idx,
+                                                            direction dir,
+                                                            bool check_collision = true) const -> std::optional<puzzle>;
+    [[nodiscard]] static auto sorted_replace(std::array<uint16_t, MAX_BLOCKS> blocks,
+                                             uint8_t idx,
+                                             uint16_t new_val) -> std::array<uint16_t, MAX_BLOCKS>;
+    [[nodiscard]] auto blocks_bitmap() const -> uint64_t;
+    [[nodiscard]] auto blocks_bitmap_h() const -> uint64_t;
+    [[nodiscard]] auto blocks_bitmap_v() const -> uint64_t;
+    static INLINE inline auto bitmap_clear_bit(uint64_t& bitmap, uint8_t w, uint8_t x, uint8_t y) -> void;
+    static INLINE inline auto bitmap_set_bit(uint64_t& bitmap, uint8_t w, uint8_t x, uint8_t y) -> void;
+    [[nodiscard]] static INLINE inline auto bitmap_get_bit(uint64_t bitmap, uint8_t w, uint8_t x, uint8_t y) -> bool;
+    INLINE inline auto bitmap_clear_block(uint64_t& bitmap, block b) const -> void;
+    INLINE inline auto bitmap_set_block(uint64_t& bitmap, block b) const -> void;
+    [[nodiscard]] INLINE inline auto bitmap_is_empty(uint64_t bitmap) const -> bool;
+    [[nodiscard]] INLINE inline auto bitmap_is_full(uint64_t bitmap) const -> bool;
+
+    /**
+     * Checks if b would collide with any block on the board.
+     *
+     * @param bitmap Board occupancy map
+     * @param b Hypothetical block to check collision with
+     * @return True if b would collide with any other block on the board
+     */
+    [[nodiscard]] INLINE inline auto bitmap_check_collision(uint64_t bitmap, block b) const -> bool;
+
+    /**
+     * Checks if b would collide with any block on the board after moving in direction dir.
+     *
+     * @param bitmap Board occupancy map
+     * @param b Existing block to check collision with
+     * @param dir Direction in which the block should be moved
+     * @return True if b would collide with any other block on the board after moving in direction dir
+     */
+    [[nodiscard]] INLINE inline auto bitmap_check_collision(uint64_t bitmap, block b, direction dir) const -> bool;

    template <typename F>
-    auto for_each_adjacent(F&& callback) const -> void
+    // ReSharper disable once CppRedundantInlineSpecifier
+    INLINE inline auto for_each_adjacent(F&& callback) const -> void
    {
        const uint64_t bitmap = blocks_bitmap();
        const bool r = get_restricted();
@ -427,9 +485,444 @@ public:

    [[nodiscard]] auto explore_state_space() const
        -> std::pair<std::vector<puzzle>, std::vector<std::pair<size_t, size_t>>>;
+
+    // Determines to which cluster a puzzle belongs. Clusters are identified by the
+    // state with the numerically smallest binary representation.
+    [[nodiscard]] auto get_cluster_id_and_solution() const -> std::pair<puzzle, bool>;
+
+    [[nodiscard]] auto bitmap_find_first_empty(uint64_t bitmap, int& x, int& y) const -> bool;
+
+    static auto generate_block_sequences(
+        const boost::unordered_flat_set<block, block_hasher2, block_equal2>& permitted_blocks,
+        block target_block,
+        size_t max_blocks,
+        std::vector<block>& current_sequence,
+        int current_area,
+        int board_area,
+        const std::function<void(const std::vector<block>&)>& callback) -> void;
+
+    static auto place_block_sequence(const puzzle& p,
+                                     const uint64_t& bitmap,
+                                     const std::tuple<uint8_t, uint8_t, uint8_t, uint8_t, bool, bool>& p_repr,
+                                     const std::vector<block>& sequence,
+                                     block target_block,
+                                     const std::tuple<uint8_t, uint8_t, uint8_t, uint8_t>& target_block_pos_range,
+                                     bool has_target,
+                                     size_t index,
+                                     const std::function<void(const puzzle&)>& callback) -> void;
+
+    [[nodiscard]] auto explore_puzzle_space(
+        const boost::unordered_flat_set<block, block_hasher2, block_equal2>& permitted_blocks,
+        block target_block,
+        const std::tuple<uint8_t, uint8_t, uint8_t, uint8_t>& target_block_pos_range,
+        size_t max_blocks,
+        std::optional<BS::thread_pool<>* const> thread_pool = std::nullopt) const -> boost::unordered_flat_set<
+        puzzle, puzzle_hasher>;
 };

-// Hash functions for sets and maps
+// Inline functions definitions
+#ifndef REGION_INLINE_DEFS
+
+inline auto puzzle::block::set_x(const uint8_t x) const -> block
+{
+    #ifdef RUNTIME_CHECKS
+    if (x > 7) {
+        throw std::invalid_argument("Block x-position out of bounds");
+    }
+    #endif
+
+    block b = *this;
+    set_bits(b.repr, X_S, X_E, x);
+    return b;
+}
+
+inline auto puzzle::block::set_y(const uint8_t y) const -> block
+{
+    #ifdef RUNTIME_CHECKS
+    if (y > 7) {
+        throw std::invalid_argument("Block y-position out of bounds");
+    }
+    #endif
+
+    block b = *this;
+    set_bits(b.repr, Y_S, Y_E, y);
+    return b;
+}
+
+inline auto puzzle::block::set_width(const uint8_t width) const -> block
+{
+    #ifdef RUNTIME_CHECKS
+    if (width - 1 > 7) {
+        throw std::invalid_argument("Block width out of bounds");
+    }
+    #endif
+
+    block b = *this;
+    set_bits(b.repr, WIDTH_S, WIDTH_E, width - 1u);
+    return b;
+}
+
+inline auto puzzle::block::set_height(const uint8_t height) const -> block
+{
+    #ifdef RUNTIME_CHECKS
+    if (height - 1 > 7) {
+        throw std::invalid_argument("Block height out of bounds");
+    }
+    #endif
+
+    block b = *this;
+    set_bits(b.repr, HEIGHT_S, HEIGHT_E, height - 1u);
+    return b;
+}
+
+inline auto puzzle::block::set_target(const bool target) const -> block
+{
+    block b = *this;
+    set_bits(b.repr, TARGET_S, TARGET_E, target);
+    return b;
+}
+
+inline auto puzzle::block::set_immovable(const bool immovable) const -> block
+{
+    block b = *this;
+    set_bits(b.repr, IMMOVABLE_S, IMMOVABLE_E, immovable);
+    return b;
+}
+
+inline auto puzzle::block::get_x() const -> uint8_t
+{
+    return get_bits(repr, X_S, X_E);
+}
+
+inline auto puzzle::block::get_y() const -> uint8_t
+{
+    return get_bits(repr, Y_S, Y_E);
+}
+
+inline auto puzzle::block::get_width() const -> uint8_t
+{
+    return get_bits(repr, WIDTH_S, WIDTH_E) + 1u;
+}
+
+inline auto puzzle::block::get_height() const -> uint8_t
+{
+    return get_bits(repr, HEIGHT_S, HEIGHT_E) + 1u;
+}
+
+inline auto puzzle::block::get_target() const -> bool
+{
+    return get_bits(repr, TARGET_S, TARGET_E);
+}
+
+inline auto puzzle::block::get_immovable() const -> bool
+{
+    return get_bits(repr, IMMOVABLE_S, IMMOVABLE_E);
+}
+
+inline auto puzzle::set_restricted(const bool restricted) const -> puzzle
+{
+    uint16_t meta = repr.cooked.meta;
+    set_bits(meta, RESTRICTED_S, RESTRICTED_E, restricted);
+    return puzzle(meta, repr.cooked.blocks);
+}
+
+inline auto puzzle::set_width(const uint8_t width) const -> puzzle
+{
+    #ifdef RUNTIME_CHECKS
+    if (width - 1 > MAX_WIDTH) {
+        throw "Board width out of bounds";
+    }
+    #endif
+
+    uint16_t meta = repr.cooked.meta;
+    set_bits(meta, WIDTH_S, WIDTH_E, width - 1u);
+    return puzzle(meta, repr.cooked.blocks);
+}
+
+inline auto puzzle::set_height(const uint8_t height) const -> puzzle
+{
+    #ifdef RUNTIME_CHECKS
+    if (height - 1 > MAX_HEIGHT) {
+        throw "Board height out of bounds";
+    }
+    #endif
+
+    uint16_t meta = repr.cooked.meta;
+    set_bits(meta, HEIGHT_S, HEIGHT_E, height - 1u);
+    return puzzle(meta, repr.cooked.blocks);
+}
+
+inline auto puzzle::set_goal(const bool goal) const -> puzzle
+{
+    uint16_t meta = repr.cooked.meta;
+    set_bits(meta, GOAL_S, GOAL_E, goal);
+    return puzzle(meta, repr.cooked.blocks);
+}
+
+inline auto puzzle::set_goal_x(const uint8_t target_x) const -> puzzle
+{
+    #ifdef RUNTIME_CHECKS
+    if (target_x >= MAX_WIDTH) {
+        throw "Board target x out of bounds";
+    }
+    #endif
+
+    uint16_t meta = repr.cooked.meta;
+    set_bits(meta, GOAL_X_S, GOAL_X_E, target_x);
+    return puzzle(meta, repr.cooked.blocks);
+}
+
+inline auto puzzle::set_goal_y(const uint8_t target_y) const -> puzzle
+{
+    #ifdef RUNTIME_CHECKS
+    if (target_y >= MAX_HEIGHT) {
+        throw "Board target y out of bounds";
+    }
+    #endif
+
+    uint16_t meta = repr.cooked.meta;
+    set_bits(meta, GOAL_Y_S, GOAL_Y_E, target_y);
+    return puzzle(meta, repr.cooked.blocks);
+}
+
+inline auto puzzle::get_restricted() const -> bool
+{
+    return get_bits(repr.cooked.meta, RESTRICTED_S, RESTRICTED_E);
+}
+
+inline auto puzzle::get_width() const -> uint8_t
+{
+    return get_bits(repr.cooked.meta, WIDTH_S, WIDTH_E) + 1u;
+}
+
+inline auto puzzle::get_height() const -> uint8_t
+{
+    return get_bits(repr.cooked.meta, HEIGHT_S, HEIGHT_E) + 1u;
+}
+
+inline auto puzzle::get_goal() const -> bool
+{
+    return get_bits(repr.cooked.meta, GOAL_S, GOAL_E);
+}
+
+inline auto puzzle::get_goal_x() const -> uint8_t
+{
+    return get_bits(repr.cooked.meta, GOAL_X_S, GOAL_X_E);
+}
+
+inline auto puzzle::get_goal_y() const -> uint8_t
+{
+    return get_bits(repr.cooked.meta, GOAL_Y_S, GOAL_Y_E);
+}
+
+INLINE inline auto puzzle::try_move_block_at_fast(uint64_t bitmap,
+                                                  const uint8_t block_idx,
+                                                  const direction dir,
+                                                  const bool check_collision) const -> std::optional<puzzle>
+{
+    const block b = block(repr.cooked.blocks[block_idx]);
+    const auto [bx, by, bw, bh, bt, bi] = b.unpack_repr();
+    if (bi) {
+        return std::nullopt;
+    }
+
+    const auto [w, h, gx, gy, r, g] = unpack_meta();
+    const int dirs = r ? b.principal_dirs() : nor | eas | sou | wes;
+
+    // Get target block
+    int _target_x = bx;
+    int _target_y = by;
+    switch (dir) {
+    case nor:
+        if (!(dirs & nor) || _target_y < 1) {
+            return std::nullopt;
+        }
+        --_target_y;
+        break;
+    case eas:
+        if (!(dirs & eas) || _target_x + bw >= w) {
+            return std::nullopt;
+        }
+        ++_target_x;
+        break;
+    case sou:
+        if (!(dirs & sou) || _target_y + bh >= h) {
+            return std::nullopt;
+        }
+        ++_target_y;
+        break;
+    case wes:
+        if (!(dirs & wes) || _target_x < 1) {
+            return std::nullopt;
+        }
+        --_target_x;
+        break;
+    }
+
+    // Check collisions
+    if (check_collision) {
+        bitmap_clear_block(bitmap, b);
+        if (bitmap_check_collision(bitmap, b, dir)) {
+            return std::nullopt;
+        }
+    }
+
+    // Replace block
+    const std::array<uint16_t, MAX_BLOCKS> blocks = sorted_replace(repr.cooked.blocks,
+                                                                   block_idx,
+                                                                   block::create_repr(
+                                                                       _target_x,
+                                                                       _target_y,
+                                                                       bw,
+                                                                       bh,
+                                                                       bt));
+
+    // This constructor doesn't sort
+    return puzzle(std::make_tuple(w, h, gx, gy, r, g), blocks);
+}
+
+INLINE inline auto puzzle::bitmap_clear_bit(uint64_t& bitmap, const uint8_t w, const uint8_t x, const uint8_t y) -> void
+{
+    set_bits(bitmap, y * w + x, y * w + x, 0u);
+}
+
+INLINE inline auto puzzle::bitmap_set_bit(uint64_t& bitmap, const uint8_t w, const uint8_t x, const uint8_t y) -> void
+{
+    set_bits(bitmap, y * w + x, y * w + x, 1u);
+}
+
+INLINE inline auto puzzle::bitmap_get_bit(const uint64_t bitmap,
+                                          const uint8_t w,
+                                          const uint8_t x,
+                                          const uint8_t y) -> bool
+{
+    return get_bits(bitmap, y * w + x, y * w + x);
+}
+
+INLINE inline auto puzzle::bitmap_clear_block(uint64_t& bitmap, const block b) const -> void
+{
+    const auto [x, y, w, h, t, i] = b.unpack_repr();
+    const uint8_t width = get_width();
+
+    for (int dy = 0; dy < h; ++dy) {
+        for (int dx = 0; dx < w; ++dx) {
+            bitmap_clear_bit(bitmap, width, x + dx, y + dy);
+        }
+    }
+}
+
+INLINE inline auto puzzle::bitmap_set_block(uint64_t& bitmap, const block b) const -> void
+{
+    const auto [x, y, w, h, t, i] = b.unpack_repr();
+    const uint8_t width = get_width();
+
+    for (int dy = 0; dy < h; ++dy) {
+        for (int dx = 0; dx < w; ++dx) {
+            bitmap_set_bit(bitmap, width, x + dx, y + dy);
+        }
+    }
+}
+
+INLINE inline auto puzzle::bitmap_is_empty(const uint64_t bitmap) const -> bool
+{
+    const uint8_t shift = 64 - get_width() * get_height();
+    return bitmap << shift == 0;
+}
+
+INLINE inline auto puzzle::bitmap_is_full(const uint64_t bitmap) const -> bool
+{
+    const uint8_t shift = 64 - get_width() * get_height();
+    return ((bitmap << shift) >> shift) == ((static_cast<uint64_t>(-1) << shift) >> shift);
+}
+
+INLINE inline auto puzzle::bitmap_check_collision(const uint64_t bitmap, const block b) const -> bool
+{
+    const auto [x, y, w, h, t, i] = b.unpack_repr();
+    const uint8_t width = get_width();
+
+    for (int dy = 0; dy < h; ++dy) {
+        for (int dx = 0; dx < w; ++dx) {
+            if (bitmap_get_bit(bitmap, width, x + dx, y + dy)) {
+                return true; // collision
+            }
+        }
+    }
+
+    return false;
+}
+
+INLINE inline auto puzzle::bitmap_check_collision(const uint64_t bitmap,
+                                                  const block b,
+                                                  const direction dir) const -> bool
+{
+    const auto [x, y, w, h, t, i] = b.unpack_repr();
+    const uint8_t width = get_width();
+
+    switch (dir) {
+    case nor: // Check the row above: (x...x+w-1, y-1)
+        for (int dx = 0; dx < w; ++dx) {
+            if (bitmap_get_bit(bitmap, width, x + dx, y - 1)) {
+                return true;
+            }
+        }
+        break;
+    case sou: // Check the row below: (x...x+w-1, y+h)
+        for (int dx = 0; dx < w; ++dx) {
+            if (bitmap_get_bit(bitmap, width, x + dx, y + h)) {
+                return true;
+            }
+        }
+        break;
+    case wes: // Check the column left: (x-1, y...y+h-1)
+        for (int dy = 0; dy < h; ++dy) {
+            if (bitmap_get_bit(bitmap, width, x - 1, y + dy)) {
+                return true;
+            }
+        }
+        break;
+    case eas: // Check the column right: (x+w, y...y+h-1)
+        for (int dy = 0; dy < h; ++dy) {
+            if (bitmap_get_bit(bitmap, width, x + w, y + dy)) {
+                return true;
+            }
+        }
+        break;
+    }
+    return false;
+}
+
+#endif
+
+// Hash functions for sets and maps.
+// Declared after puzzle class to use puzzle::hash_combine
+#ifndef REGION_HASHERS
+
+struct block_hasher
+{
+    auto operator()(const puzzle::block& b) const noexcept -> size_t
+    {
+        return b.hash();
+    }
+};
+
+struct block_hasher2
+{
+    auto operator()(const puzzle::block& b) const noexcept -> size_t
+    {
+        return b.position_independent_hash();
+    }
+};
+
+struct block_equal2
+{
+    auto operator()(const puzzle::block& a, const puzzle::block& b) const noexcept -> bool
+    {
+        const auto [ax, ay, aw, ah, at, ai] = a.unpack_repr();
+        const auto [bx, by, bw, bh, bt, bi] = b.unpack_repr();
+
+        return aw == bw && ah == bh && at == bt && ai == bi;
+    }
+};

 struct puzzle_hasher
 {
@ -453,7 +946,7 @@ struct link_hasher
    }
 };

-struct link_equal_to
+struct link_equal
 {
    auto operator()(const std::pair<puzzle, puzzle>& a, const std::pair<puzzle, puzzle>& b) const noexcept -> bool
    {
@ -462,9 +955,9 @@ struct link_equal_to
 };

 template <typename T, typename... Rest>
-auto puzzle::hash_combine(std::size_t& seed, const T& v, const Rest&... rest) -> void
+auto puzzle::hash_combine(size_t& seed, const T& v, const Rest&... rest) -> void
 {
-    auto h = []<typename HashedType>(const HashedType& val) -> std::size_t
+    auto hasher = []<typename HashedType>(const HashedType& val) -> std::size_t
    {
        if constexpr (std::is_same_v<std::decay_t<HashedType>, puzzle>) {
            return puzzle_hasher{}(val);
@ -475,8 +968,10 @@ auto puzzle::hash_combine(std::size_t& seed, const T& v, const Rest&... rest) ->
        }
    };

-    seed ^= h(v) + 0x9e3779b9 + (seed << 6) + (seed >> 2);
+    seed ^= hasher(v) + 0x9e3779b9 + (seed << 6) + (seed >> 2);
    (hash_combine(seed, rest), ...);
 }

+#endif
+
 #endif
--- a/include/renderer.hpp
+++ b/include/renderer.hpp
@ -1,8 +1,8 @@
 #ifndef RENDERER_HPP_
 #define RENDERER_HPP_

-#include "orbit_camera.hpp"
 #include "config.hpp"
+#include "orbit_camera.hpp"
 #include "input_handler.hpp"
 #include "state_manager.hpp"
 #include "user_interface.hpp"
@ -18,7 +18,7 @@ private:
    user_interface& gui;

    const orbit_camera& camera;
-    RenderTexture render_target =
+    RenderTexture graph_target =
        LoadRenderTexture(GetScreenWidth() / 2, GetScreenHeight() - MENU_HEIGHT);

    // TODO: Those should be moved to the user_interface.h
@ -81,7 +81,7 @@ public:

    ~renderer()
    {
-        UnloadRenderTexture(render_target);
+        UnloadRenderTexture(graph_target);
        UnloadRenderTexture(klotski_target);
        UnloadRenderTexture(menu_target);

--- a/include/state_manager.hpp
+++ b/include/state_manager.hpp
@ -2,10 +2,10 @@
 #define STATE_MANAGER_HPP_

 #include "graph_distances.hpp"
+#include "load_save.hpp"
 #include "threaded_physics.hpp"
 #include "puzzle.hpp"

-#include <stack>
 #include <boost/unordered/unordered_flat_map.hpp>
 #include <boost/unordered/unordered_flat_set.hpp>

@ -43,10 +43,9 @@ private:

 public:
    state_manager(threaded_physics& _physics, const std::string& _preset_file)
-        : physics(_physics)
+        : physics(_physics), preset_file(_preset_file)
    {
-        parse_preset_file(_preset_file);
-        load_preset(0);
+        reload_preset_file();
    }

    state_manager(const state_manager& copy) = delete;
@ -95,15 +94,13 @@ private:

 public:
    // Presets
-
-    auto parse_preset_file(const std::string& _preset_file) -> bool;
-    auto append_preset_file(const std::string& preset_name) -> bool;
+    auto save_current_to_preset_file(const std::string& preset_comment) -> void;
+    auto reload_preset_file() -> void;
    auto load_preset(size_t preset) -> void;
    auto load_previous_preset() -> void;
    auto load_next_preset() -> void;

    // Update current_state
-
    auto update_current_state(const puzzle& p) -> void;
    auto edit_starting_state(const puzzle& p) -> void;
    auto goto_starting_state() -> void;
@ -113,7 +110,6 @@ public:
    auto goto_closest_target_state() -> void;

    // Update graph
-
    auto populate_graph() -> void;
    auto clear_graph_and_add_current(const puzzle& p) -> void;
    auto clear_graph_and_add_current() -> void;
@ -122,7 +118,6 @@ public:
    auto populate_winning_path() -> void;

    // Index mapping
-
    [[nodiscard]] auto get_index(const puzzle& state) const -> size_t;
    [[nodiscard]] auto get_current_index() const -> size_t;
    [[nodiscard]] auto get_starting_index() const -> size_t;
--- a/include/threaded_physics.hpp
+++ b/include/threaded_physics.hpp
@ -1,6 +1,8 @@
 #ifndef PHYSICS_HPP_
 #define PHYSICS_HPP_

+#include "config.hpp"
+
 #include <atomic>
 #include <condition_variable>
 #include <mutex>
@ -11,13 +13,10 @@
 #include <variant>
 #include <vector>

-#ifdef TRACY
-#include <tracy/Tracy.hpp>
-#endif
-
 class threaded_physics
 {
-    struct add_mass {};
+    struct add_mass
+    {};

    struct add_spring
    {
@ -25,24 +24,25 @@ class threaded_physics
        size_t b;
    };

-    struct clear_graph {};
+    struct clear_graph
+    {};

    using command = std::variant<add_mass, add_spring, clear_graph>;

    struct physics_state
    {
-        #ifdef TRACY
+#ifdef TRACY
        TracyLockable(std::mutex, command_mtx);
-        #else
+#else
        std::mutex command_mtx;
-        #endif
+#endif
        std::queue<command> pending_commands;

-        #ifdef TRACY
+#ifdef TRACY
        TracyLockable(std::mutex, data_mtx);
-        #else
+#else
        std::mutex data_mtx;
-        #endif
+#endif
        std::condition_variable_any data_ready_cnd;
        std::condition_variable_any data_consumed_cnd;
        Vector3 mass_center = Vector3Zero();
@ -57,14 +57,17 @@ class threaded_physics
    };

 private:
+    std::optional<BS::thread_pool<>* const> thread_pool;
    std::thread physics;

 public:
    physics_state state;

 public:
-    threaded_physics()
-        : physics(physics_thread, std::ref(state)) {}
+    explicit threaded_physics(
+        const std::optional<BS::thread_pool<>* const> _thread_pool = std::nullopt)
+        : thread_pool(_thread_pool), physics(physics_thread, std::ref(state), std::ref(thread_pool))
+    {}

    threaded_physics(const threaded_physics& copy) = delete;
    auto operator=(const threaded_physics& copy) -> threaded_physics& = delete;
@ -80,16 +83,19 @@ public:
    }

 private:
-    static auto physics_thread(physics_state& state) -> void;
+#ifdef ASYNC_OCTREE
+    static auto set_octree_pool_thread_name(size_t idx) -> void;
+#endif
+
+    static auto physics_thread(physics_state& state,
+                               std::optional<BS::thread_pool<>* const> thread_pool) -> void;

 public:
-    auto add_mass_cmd() -> void;
-
-    auto add_spring_cmd(size_t a, size_t b) -> void;
-
    auto clear_cmd() -> void;
-
-    auto add_mass_springs_cmd(size_t num_masses, const std::vector<std::pair<size_t, size_t>>& springs) -> void;
+    auto add_mass_cmd() -> void;
+    auto add_spring_cmd(size_t a, size_t b) -> void;
+    auto add_mass_springs_cmd(size_t num_masses,
+                              const std::vector<std::pair<size_t, size_t>>& springs) -> void;
 };

-#endif
+#endif
--- a/include/user_interface.hpp
+++ b/include/user_interface.hpp
@ -1,5 +1,5 @@
-#ifndef GUI_HPP_
-#define GUI_HPP_
+#ifndef USER_INTERFACE_HPP_
+#define USER_INTERFACE_HPP_

 #include "orbit_camera.hpp"
 #include "config.hpp"
@ -104,7 +104,7 @@ private:
    bool ok_message = false;
    bool yes_no_message = false;
    bool save_window = false;
-    std::array<char, 256> preset_name = {};
+    std::array<char, 256> preset_comment = {};
    bool help_window = false;

 public:
--- a/include/util.hpp
+++ b/include/util.hpp
@ -4,81 +4,8 @@
 #include <iostream>
 #include <raylib.h>

-// Bit shifting + masking
-
-template <class T>
-    requires std::unsigned_integral<T>
-auto create_mask(const uint8_t first, const uint8_t last) -> T
-{
-    // If the mask width is equal the type width return all 1s instead of shifting
-    // as shifting by type-width is undefined behavior.
-    if (static_cast<size_t>(last - first + 1) >= sizeof(T) * 8) {
-        return ~T{0};
-    }
-
-    // Example: first=4, last=7, 7-4+1=4
-    //           1 << 4 = 0b00010000
-    //          32  - 1 = 0b00001111
-    //          31 << 4 = 0b11110000
-    // Subtracting 1 generates a consecutive mask.
-    return ((T{1} << (last - first + 1)) - 1) << first;
-}
-
-template <class T>
-    requires std::unsigned_integral<T>
-auto clear_bits(T& bits, const uint8_t first, const uint8_t last) -> void
-{
-    const T mask = create_mask<T>(first, last);
-
-    bits = bits & ~mask;
-}
-
-template <class T, class U>
-    requires std::unsigned_integral<T> && std::unsigned_integral<U>
-auto set_bits(T& bits, const uint8_t first, const uint8_t last, const U value) -> void
-{
-    const T mask = create_mask<T>(first, last);
-
-    // Example: first=4, last=6, value=0b1110,   bits = 0b 01111110
-    //          mask                                  = 0b 01110000
-    //          bits & ~mask                          = 0b 00001110
-    //          value << 4                            = 0b 11100000
-    //          (value << 4) & mask                   = 0b 01100000
-    //          (bits & ~mask) | (value << 4) & mask  = 0b 01101110
-    //                                 Insert position:     ^^^
-    // First clear the bits, then | with the value positioned at the insertion point.
-    // The value may be larger than [first, last], extra bits are ignored.
-    bits = (bits & ~mask) | ((static_cast<T>(value) << first) & mask);
-}
-
-template <class T>
-    requires std::unsigned_integral<T>
-auto get_bits(const T bits, const uint8_t first, const uint8_t last) -> T
-{
-    const T mask = create_mask<T>(first, last);
-
-    // We can >> without sign extension because T is unsigned_integral
-    return (bits & mask) >> first;
-}
-
-// std::variant visitor
-
-// https://en.cppreference.com/w/cpp/utility/variant/visit
-template <class... Ts>
-struct overloads : Ts...
-{
-    using Ts::operator()...;
-};
-
-// Enums
-
-enum direction
-{
-    nor = 1 << 0,
-    eas = 1 << 1,
-    sou = 1 << 2,
-    wes = 1 << 3,
-};
+#define INLINE __attribute__((always_inline))
+#define PACKED __attribute__((packed))

 enum ctrl
 {
@ -115,6 +42,96 @@ enum bg
    bg_white = 47
 };

+inline auto ansi_bold_fg(const fg color) -> std::string
+{
+    return std::format("\033[{};{}m", static_cast<int>(bold_bright), static_cast<int>(color));
+}
+
+inline auto ansi_reset() -> std::string
+{
+    return std::format("\033[{}m", static_cast<int>(reset));
+}
+
+// Bit shifting + masking
+
+template <class T>
+    requires std::unsigned_integral<T>
+// ReSharper disable once CppRedundantInlineSpecifier
+INLINE inline auto create_mask(const uint8_t first, const uint8_t last) -> T
+{
+    // If the mask width is equal the type width return all 1s instead of shifting
+    // as shifting by type-width is undefined behavior.
+    if (static_cast<size_t>(last - first + 1) >= sizeof(T) * 8) {
+        return ~T{0};
+    }
+
+    // Example: first=4, last=7, 7-4+1=4
+    //           1 << 4 = 0b00010000
+    //          32  - 1 = 0b00001111
+    //          31 << 4 = 0b11110000
+    // Subtracting 1 generates a consecutive mask.
+    return ((T{1} << (last - first + 1)) - 1) << first;
+}
+
+template <class T>
+    requires std::unsigned_integral<T>
+// ReSharper disable once CppRedundantInlineSpecifier
+INLINE inline auto clear_bits(T& bits, const uint8_t first, const uint8_t last) -> void
+{
+    const T mask = create_mask<T>(first, last);
+
+    bits = bits & ~mask;
+}
+
+template <class T, class U>
+    requires std::unsigned_integral<T> && std::unsigned_integral<U>
+// ReSharper disable once CppRedundantInlineSpecifier
+INLINE inline auto set_bits(T& bits, const uint8_t first, const uint8_t last, const U value) -> void
+{
+    const T mask = create_mask<T>(first, last);
+
+    // Example: first=4, last=6, value=0b1110,   bits = 0b 01111110
+    //          mask                                  = 0b 01110000
+    //          bits & ~mask                          = 0b 00001110
+    //          value << 4                            = 0b 11100000
+    //          (value << 4) & mask                   = 0b 01100000
+    //          (bits & ~mask) | (value << 4) & mask  = 0b 01101110
+    //                                 Insert position:     ^^^
+    // First clear the bits, then | with the value positioned at the insertion point.
+    // The value may be larger than [first, last], extra bits are ignored.
+    bits = (bits & ~mask) | ((static_cast<T>(value) << first) & mask);
+}
+
+template <class T>
+    requires std::unsigned_integral<T>
+// ReSharper disable once CppRedundantInlineSpecifier
+INLINE inline auto get_bits(const T bits, const uint8_t first, const uint8_t last) -> T
+{
+    const T mask = create_mask<T>(first, last);
+
+    // We can >> without sign extension because T is unsigned_integral
+    return (bits & mask) >> first;
+}
+
+// std::variant visitor
+
+// https://en.cppreference.com/w/cpp/utility/variant/visit
+template <class... Ts>
+struct overloads : Ts...
+{
+    using Ts::operator()...;
+};
+
+// Enums
+
+enum direction
+{
+    nor = 1 << 0,
+    eas = 1 << 1,
+    sou = 1 << 2,
+    wes = 1 << 3,
+};
+
 // Output

 inline auto operator<<(std::ostream& os, const Vector2& v) -> std::ostream&
@ -129,17 +146,14 @@ inline auto operator<<(std::ostream& os, const Vector3& v) -> std::ostream&
    return os;
 }

-inline auto ansi_bold_fg(const fg color) -> std::string
-{
-    return std::format("\033[1;{}m", static_cast<int>(color));
-}
-
-inline auto ansi_reset() -> std::string
-{
-    return "\033[0m";
-}
-
 // std::println doesn't work with mingw
+template <typename... Args>
+auto traceln(std::format_string<Args...> fmt, Args&&... args) -> void
+{
+    std::cout << std::format("[{}TRACE{}]: ", ansi_bold_fg(fg_cyan), ansi_reset()) << std::format(
+        fmt, std::forward<Args>(args)...) << std::endl;
+}
+
 template <typename... Args>
 auto infoln(std::format_string<Args...> fmt, Args&&... args) -> void
 {
@ -161,4 +175,19 @@ auto errln(std::format_string<Args...> fmt, Args&&... args) -> void
        fmt, std::forward<Args>(args)...) << std::endl;
 }

+inline auto print_bitmap(const uint64_t bitmap, const uint8_t w, const uint8_t h, const std::string& title) -> void
+{
+    traceln("{}:", title);
+    traceln("{}", std::string(2 * w - 1, '='));
+    for (size_t y = 0; y < w; ++y) {
+        std::cout << "         ";
+        for (size_t x = 0; x < h; ++x) {
+            std::cout << static_cast<int>(get_bits(bitmap, y * w + x, y * h + x)) << " ";
+        }
+        std::cout << "\n";
+    }
+    std::cout << std::flush;
+    traceln("{}", std::string(2 * w - 1, '='));
+}
+
 #endif
--- a/src/graph_distances.cpp
+++ b/src/graph_distances.cpp
@ -2,10 +2,6 @@

 #include <queue>

-#ifdef TRACY
-    #include <tracy/Tracy.hpp>
-#endif
-
 auto graph_distances::clear() -> void
 {
    distances.clear();
--- a/src/input_handler.cpp
+++ b/src/input_handler.cpp
@ -278,7 +278,7 @@ auto input_handler::toggle_camera_projection() const -> void
 auto input_handler::move_block_nor() -> void
 {
    const puzzle& current = state.get_current_state();
-    const std::optional<puzzle>& next = current.try_move_block_at_fast(sel_x, sel_y, nor);
+    const std::optional<puzzle>& next = current.try_move_block_at(sel_x, sel_y, nor);
    if (!next) {
        return;
    }
@ -290,7 +290,7 @@ auto input_handler::move_block_nor() -> void
 auto input_handler::move_block_wes() -> void
 {
    const puzzle& current = state.get_current_state();
-    const std::optional<puzzle>& next = current.try_move_block_at_fast(sel_x, sel_y, wes);
+    const std::optional<puzzle>& next = current.try_move_block_at(sel_x, sel_y, wes);
    if (!next) {
        return;
    }
@ -302,7 +302,7 @@ auto input_handler::move_block_wes() -> void
 auto input_handler::move_block_sou() -> void
 {
    const puzzle& current = state.get_current_state();
-    const std::optional<puzzle>& next = current.try_move_block_at_fast(sel_x, sel_y, sou);
+    const std::optional<puzzle>& next = current.try_move_block_at(sel_x, sel_y, sou);
    if (!next) {
        return;
    }
@ -314,7 +314,7 @@ auto input_handler::move_block_sou() -> void
 auto input_handler::move_block_eas() -> void
 {
    const puzzle& current = state.get_current_state();
-    const std::optional<puzzle>& next = current.try_move_block_at_fast(sel_x, sel_y, eas);
+    const std::optional<puzzle>& next = current.try_move_block_at(sel_x, sel_y, eas);
    if (!next) {
        return;
    }
--- a/src/load_save.cpp
+++ b/src/load_save.cpp
@ -0,0 +1,79 @@
+#include "load_save.hpp"
+
+#include <fstream>
+
+auto parse_preset_file(const std::string& preset_file) -> std::pair<std::vector<puzzle>, std::vector<std::string>>
+{
+    std::fstream file(preset_file, std::ios::in);
+    if (!file) {
+        infoln("Preset file \"{}\" couldn't be opened.", preset_file);
+        return {};
+    }
+
+    std::string line;
+    std::vector<std::string> comment_lines;
+    std::vector<std::string> preset_lines;
+    while (std::getline(file, line)) {
+        if (line.starts_with("S")) {
+            preset_lines.push_back(line);
+        } else if (line.starts_with("#")) {
+            comment_lines.push_back(line);
+        }
+    }
+
+    if (preset_lines.empty() || comment_lines.size() != preset_lines.size()) {
+        infoln("Preset file \"{}\" couldn't be opened.", preset_file);
+        return {};
+    }
+
+    std::vector<puzzle> preset_states;
+    for (const auto& preset : preset_lines) {
+        // Each char is a bit
+        const puzzle& p = puzzle(preset);
+
+        if (const std::optional<std::string>& reason = p.try_get_invalid_reason()) {
+            infoln("Preset file \"{}\" contained invalid presets: {}", preset_file, *reason);
+            return {};
+        }
+        preset_states.emplace_back(p);
+    }
+
+    infoln("Loaded {} presets from \"{}\".", preset_lines.size(), preset_file);
+
+    return {preset_states, comment_lines};
+}
+
+auto append_preset_file(const std::string& preset_file, const std::string& preset_name, const puzzle& p) -> bool
+{
+    infoln(R"(Saving preset "{}" to "{}")", preset_name, preset_file);
+
+    if (p.try_get_invalid_reason()) {
+        return false;
+    }
+
+    std::fstream file(preset_file, std::ios_base::app | std::ios_base::out);
+    if (!file) {
+        infoln("Preset file \"{}\" couldn't be opened.", preset_file);
+        return false;
+    }
+
+    file << "\n# " << preset_name << "\n" << p.string_repr() << std::flush;
+
+    return true;
+}
+
+auto append_preset_file_quiet(const std::string& preset_file, const std::string& preset_name, const puzzle& p, const bool validate) -> bool
+{
+    if (validate && p.try_get_invalid_reason()) {
+        return false;
+    }
+
+    std::fstream file(preset_file, std::ios_base::app | std::ios_base::out);
+    if (!file) {
+        return false;
+    }
+
+    file << "\n# " << preset_name << "\n" << p.string_repr() << std::flush;
+
+    return true;
+}
--- a/src/main.cpp
+++ b/src/main.cpp
@ -4,75 +4,60 @@

 #include "config.hpp"
 #include "input_handler.hpp"
-#include "mass_spring_system.hpp"
 #include "threaded_physics.hpp"
 #include "renderer.hpp"
 #include "state_manager.hpp"
 #include "user_interface.hpp"

-#ifdef TRACY
-#include <tracy/Tracy.hpp>
+#include <filesystem>
+
+#ifndef WIN32
+#include <boost/program_options.hpp>
+namespace po = boost::program_options;
 #endif

-// TODO: Add some popups (my split between input.cpp/gui.cpp makes this ugly)
-//       - Clear graph: Notify that this will clear the visited states and move
-//         history
-//       - Reset state: Notify that this will reset the move count

-// TODO: Reduce memory usage
-//       - The memory model of the puzzle board is terrible (bitboards?)
-
-// TODO: Improve solver
-//       - Move discovery is terrible
-//         - Instead of trying each direction for each block, determine the
-//           possible moves more efficiently (requires a different memory model)
-//       - Implement state discovery/enumeration
-//         - Find all possible initial board states (single one for each
-//           possible statespace). Currently wer're just finding all states
-//           given the initial state
-//         - Would allow to generate random puzzles with a certain move count
-
-// TODO: Move selection accordingly when undoing moves (need to diff two states
-//       and get the moved blocks)
+// TODO: Implement state discovery/enumeration
+//       - Find all possible initial board states (single one for each possible statespace).
+//         Currently wer're just finding all states given the initial state
+//       - Would allow to generate random puzzles with a certain move count
+// TODO: Export cluster to graphviz
+// TODO: Fix naming:
+//       - Target: The block that has to leave the board to win
+//       - Goal: The opening in the board for the target
+//       - Puzzle (not board or state): A puzzle configuration (width, height, goal_x, goal_y, restricted, goal)
+//       - Block: A puzzle block (x, y, width, height, target, immovable)
+//       - Puzzle State: A specific puzzle state (width, height, goal_x, goal_y, restricted, goal, blocks)
+//       - Cluster: A graph of puzzle states connected by moves, generated from a specific Puzzle State
+//       - Puzzle Space: A number of Clusters generated from a generic Puzzle
+// TODO: Add state space generation time to debug overlay
+// TODO: Move selection accordingly when undoing moves (need to diff two states and get the moved blocks)

 // TODO: Click states in the graph to display them in the board

-// NOTE: Tracy uses a huge amount of memory. For longer testing disable Tracy.
-
-// For profiling explore_state_space
-auto main2(int argc, char* argv[]) -> int
+#ifdef THREADPOOL
+auto set_pool_thread_name(size_t idx) -> void
 {
-    const puzzle p = puzzle(
-        "S:[4x5] G:[1,3] M:[F] B:[{_ 2X2 _ _} {1x1 _ _ 1x1} {1x2 2x1 _ 1x2} {_ 2x1 _ _} {1x1 2x1 _ 1x1}]");
-
-    for (int i = 0; i < 50; ++i) {
-        auto space = p.explore_state_space();
-    }
-
-    return 0;
+    BS::this_thread::set_os_thread_name(std::format("worker-{}", idx));
 }

-auto main(int argc, char* argv[]) -> int
+BS::thread_pool<> threads(std::thread::hardware_concurrency() - 2, set_pool_thread_name);
+constexpr std::optional<BS::thread_pool<>* const> thread_pool = &threads;
+#else
+constexpr std::optional<BS::thread_pool<>* const> thread_pool = std::nullopt;
+#endif
+
+std::string preset_file;
+std::string output_file;
+int max_blocks = 5;
+int board_width = 6;
+int board_height = 6;
+int goal_x = 4;
+int goal_y = 2;
+bool restricted = true;
+
+auto ui_mode() -> int
 {
-    std::string preset_file;
-    if (argc != 2) {
-        preset_file = "default.puzzle";
-    } else {
-        preset_file = argv[1];
-    }
-
-    #ifdef BACKWARD
-    infoln("Backward stack-traces enabled.");
-    #else
-    infoln("Backward stack-traces disabled.");
-    #endif
-
-    #ifdef TRACY
-    infoln("Tracy adapter enabled.");
-    #else
-    infoln("Tracy adapter disabled.");
-    #endif
-
    // RayLib window setup
    SetTraceLogLevel(LOG_ERROR);
    SetConfigFlags(FLAG_VSYNC_HINT);
@ -82,7 +67,7 @@ auto main(int argc, char* argv[]) -> int
    InitWindow(INITIAL_WIDTH * 2, INITIAL_HEIGHT + MENU_HEIGHT, "MassSprings");

    // Game setup
-    threaded_physics physics;
+    threaded_physics physics(thread_pool);
    state_manager state(physics, preset_file);
    orbit_camera camera;
    input_handler input(state, camera);
@ -150,8 +135,8 @@ auto main(int argc, char* argv[]) -> int
        // Update the camera after the physics, so target lock is smooth
        size_t current_index = state.get_current_index();
        if (masses.size() > current_index) {
-            const mass_spring_system::mass& current_mass = mass_spring_system::mass(masses.at(current_index));
-            camera.update(current_mass.position, mass_center, input.camera_lock, input.camera_mass_center_lock);
+            const Vector3& current_mass = masses[current_index];
+            camera.update(current_mass, mass_center, input.camera_lock, input.camera_mass_center_lock);
        }

        // Rendering
@ -166,11 +151,199 @@ auto main(int argc, char* argv[]) -> int
        ++loop_iterations;

        #ifdef TRACY
-        FrameMark; FrameMarkEnd("MainThread");
+        FrameMark;
+        FrameMarkEnd("MainThread");
        #endif
    }

    CloseWindow();

    return 0;
-}
+}
+
+auto rush_hour_puzzle_space() -> int
+{
+    const boost::unordered_flat_set<puzzle::block, block_hasher2, block_equal2> permitted_blocks = {
+        puzzle::block(0, 0, 2, 1, false, false),
+        puzzle::block(0, 0, 3, 1, false, false),
+        puzzle::block(0, 0, 1, 2, false, false),
+        puzzle::block(0, 0, 1, 3, false, false)
+    };
+    const puzzle::block target_block = puzzle::block(0, 0, 2, 1, true, false);
+    const std::tuple<uint8_t, uint8_t, uint8_t, uint8_t> target_block_pos_range = {0, goal_y, board_width - 1, goal_y};
+
+    infoln("Exploring Rush-Hour puzzle space:");
+    infoln("- Size:       {}x{}", board_width, board_height);
+    infoln("- Goal:       {},{}", goal_x, goal_y);
+    infoln("- Restricted: {}", restricted);
+    infoln("- Max Blocks: {}", max_blocks);
+    infoln("- Target:     {}x{}", target_block.get_width(), target_block.get_height());
+    infoln("- Permitted block sizes:");
+    std::cout << "         ";
+    for (const puzzle::block b : permitted_blocks) {
+        std::cout << std::format(" {}x{},", b.get_width(), b.get_height());
+    }
+    std::cout << std::endl;
+    const std::chrono::high_resolution_clock::time_point start = std::chrono::high_resolution_clock::now();
+
+    const puzzle p = puzzle(board_width, board_height, goal_x, goal_y, restricted, true);
+    const boost::unordered_flat_set<puzzle, puzzle_hasher> result = p.explore_puzzle_space(
+        permitted_blocks,
+        target_block,
+        target_block_pos_range,
+        max_blocks,
+        thread_pool);
+
+    const std::chrono::high_resolution_clock::time_point end = std::chrono::high_resolution_clock::now();
+    infoln("Found {} different clusters. Took {}s.",
+           result.size(),
+           std::chrono::duration_cast<std::chrono::seconds>(end - start).count());
+
+    infoln("Sorting clusters...");
+    std::vector<puzzle> result_sorted{result.begin(), result.end()};
+    std::ranges::sort(result_sorted, std::ranges::greater{});
+    // for (const puzzle& _p : result_sorted) {
+    //     traceln("{}", _p.string_repr());
+    // }
+
+    size_t i = 0;
+    size_t success = 0;
+    std::filesystem::remove(output_file);
+    for (const puzzle& _p : result_sorted) {
+        if (append_preset_file_quiet(output_file, std::format("Cluster {}", i), _p, true)) {
+            ++success;
+        }
+        ++i;
+    }
+    if (success != result_sorted.size()) {
+        warnln("Saved {} of {} clusters", success, result_sorted.size());
+    } else {
+        infoln("Saved {} of {} clusters", success, result_sorted.size());
+    }
+
+    return 0;
+}
+
+enum class runmode
+{
+    USER_INTERFACE, RUSH_HOUR_PUZZLE_SPACE, EXIT,
+};
+
+auto argparse(const int argc, char* argv[]) -> runmode
+{
+#ifndef WIN32
+    po::options_description desc("Allowed options");
+    desc.add_options()                   //
+        ("help", "produce help message") //
+
+        ("presets", po::value<std::string>()->default_value("default.puzzle"), "load presets from file") //
+
+        ("output", po::value<std::string>()->default_value("clusters.puzzle"), "output file for generated clusters") //
+        ("space", po::value<std::string>()->value_name("rh|klotski"), "generate puzzle space with ruleset")          //
+        ("w", po::value<int>()->default_value(6), "board width")                                                     //
+        ("h", po::value<int>()->default_value(6), "board height")                                                    //
+        ("gx", po::value<int>()->default_value(4), "board goal horizontal position")                                 //
+        ("gy", po::value<int>()->default_value(2), "board goal vertical position")                                   //
+        ("free", "allow free block movement")                                                                        //
+        ("blocks", po::value<int>()->default_value(5), "block limit for puzzle space generation")                    //
+        ;
+
+    po::positional_options_description positional;
+    positional.add("presets", -1);
+
+    po::variables_map vm;
+    po::store(po::command_line_parser(argc, argv).options(desc).positional(positional).run(), vm);
+    po::notify(vm);
+
+    if (vm.contains("help")) {
+        std::cout << desc << std::endl;
+        return runmode::EXIT;
+    }
+
+    if (vm.contains("output")) {
+        output_file = vm["output"].as<std::string>();
+    }
+
+    if (vm.contains("w")) {
+        board_width = vm["w"].as<int>();
+        board_width = std::max(static_cast<int>(puzzle::MIN_WIDTH),
+                               std::min(board_width, static_cast<int>(puzzle::MAX_WIDTH)));
+    }
+
+    if (vm.contains("h")) {
+        board_height = vm["h"].as<int>();
+        board_height = std::max(static_cast<int>(puzzle::MIN_HEIGHT),
+                                std::min(board_height, static_cast<int>(puzzle::MAX_HEIGHT)));
+    }
+
+    if (vm.contains("gx")) {
+        goal_x = vm["gx"].as<int>();
+        goal_x = std::max(0, std::min(goal_x, static_cast<int>(puzzle::MAX_WIDTH) - 1));
+    }
+
+    if (vm.contains("gy")) {
+        goal_y = vm["gy"].as<int>();
+        goal_y = std::max(0, std::min(goal_y, static_cast<int>(puzzle::MAX_HEIGHT) - 1));
+    }
+
+    if (vm.contains("free")) {
+        restricted = false;
+    }
+
+    if (vm.contains("blocks")) {
+        max_blocks = vm["blocks"].as<int>();
+        max_blocks = std::max(1, std::min(max_blocks, static_cast<int>(puzzle::MAX_BLOCKS)));
+    }
+
+    if (vm.contains("space")) {
+        const std::string ruleset = vm["space"].as<std::string>();
+        if (ruleset == "rh") {
+            return runmode::RUSH_HOUR_PUZZLE_SPACE;
+        }
+        if (ruleset == "klotski") {
+            throw std::runtime_error("Not implemented");
+        }
+    }
+
+    if (vm.contains("presets")) {
+        preset_file = vm["presets"].as<std::string>();
+    }
+#endif
+
+    return runmode::USER_INTERFACE;
+}
+
+auto main(const int argc, char* argv[]) -> int
+{
+    #ifdef BACKWARD
+    infoln("Backward stack-traces enabled.");
+    #else
+    infoln("Backward stack-traces disabled.");
+    #endif
+
+    #ifdef TRACY
+    infoln("Tracy adapter enabled.");
+    #else
+    infoln("Tracy adapter disabled.");
+    #endif
+
+    infoln("Using background thread for physics.");
+    infoln("Using octree-barnes-hut for graph layout.");
+
+    #ifdef THREADPOOL
+    infoln("Additional thread-pool enabled ({} threads).", threads.get_thread_count());
+    #else
+    infoln("Additional thread-pool disabled.");
+    #endif
+
+    switch (argparse(argc, argv)) {
+    case runmode::USER_INTERFACE:
+        return ui_mode();
+    case runmode::RUSH_HOUR_PUZZLE_SPACE:
+        return rush_hour_puzzle_space();
+    case runmode::EXIT:
+        return 0;
+    };
+
+    return 1;
+}
--- a/src/mass_spring_system.cpp
+++ b/src/mass_spring_system.cpp
@ -2,66 +2,63 @@
 #include "config.hpp"

 #include <cfloat>
+#include <cstring>

-#ifdef TRACY
-#include <tracy/Tracy.hpp>
-#endif
-
-auto mass_spring_system::mass::clear_force() -> void
+auto mass_spring_system::calculate_spring_force(const size_t s) -> void
 {
-    force = Vector3Zero();
+    const spring _s = springs[s];
+    const Vector3 a_pos = positions[_s.a];
+    const Vector3 b_pos = positions[_s.b];
+    const Vector3 a_vel = velocities[_s.a];
+    const Vector3 b_vel = velocities[_s.b];
+
+    const Vector3 delta_pos = a_pos - b_pos;
+    const Vector3 delta_vel = a_vel - b_vel;
+
+    const float sq_len = Vector3DotProduct(delta_pos, delta_pos);
+    const float inv_len = 1.0f / sqrt(sq_len);
+    const float len = sq_len * inv_len;
+
+    const float hooke = SPRING_CONSTANT * (len - REST_LENGTH);
+    const float dampening = DAMPENING_CONSTANT * Vector3DotProduct(delta_vel, delta_pos) * inv_len;
+
+    const Vector3 a_force = Vector3Scale(delta_pos, -(hooke + dampening) * inv_len);
+    const Vector3 b_force = a_force * -1.0f;
+
+    forces[_s.a] += a_force;
+    forces[_s.b] += b_force;
 }

-auto mass_spring_system::mass::calculate_velocity(const float delta_time) -> void
+auto mass_spring_system::integrate_velocity(const size_t m, const float dt) -> void
 {
-    const Vector3 acceleration = Vector3Scale(force, 1.0 / MASS);
-    const Vector3 temp = Vector3Scale(acceleration, delta_time);
-    velocity = Vector3Add(velocity, temp);
+    const Vector3 acc = forces[m] / MASS;
+    velocities[m] += acc * dt;
 }

-auto mass_spring_system::mass::calculate_position(const float delta_time) -> void
+auto mass_spring_system::integrate_position(const size_t m, const float dt) -> void
 {
-    previous_position = position;
-
-    const Vector3 temp = Vector3Scale(velocity, delta_time);
-    position = Vector3Add(position, temp);
+    previous_positions[m] = positions[m];
+    positions[m] += velocities[m] * dt;
 }

-auto mass_spring_system::mass::verlet_update(const float delta_time) -> void
+auto mass_spring_system::verlet_update(const size_t m, const float dt) -> void
 {
-    const Vector3 acceleration = Vector3Scale(force, 1.0 / MASS);
-    const Vector3 temp_position = position;
+    const Vector3 acc = (forces[m] / MASS) * dt * dt;
+    const Vector3 pos = positions[m];

-    Vector3 displacement = Vector3Subtract(position, previous_position);
-    const Vector3 accel_term = Vector3Scale(acceleration, delta_time * delta_time);
+    Vector3 delta_pos = pos - previous_positions[m];
+    delta_pos *= 1.0 - VERLET_DAMPENING; // Minimal dampening

-    // Minimal dampening
-    displacement = Vector3Scale(displacement, 1.0 - VERLET_DAMPENING);
-
-    position = Vector3Add(Vector3Add(position, displacement), accel_term);
-    previous_position = temp_position;
-}
-
-auto mass_spring_system::spring::calculate_spring_force(mass& _a, mass& _b) -> void
-{
-    // TODO: Use a bungee force here instead of springs, since we already have global repulsion?
-    const Vector3 delta_position = Vector3Subtract(_a.position, _b.position);
-    const float current_length = Vector3Length(delta_position);
-    const Vector3 delta_velocity = Vector3Subtract(_a.velocity, _b.velocity);
-
-    const float hooke = SPRING_CONSTANT * (current_length - REST_LENGTH);
-    const float dampening = DAMPENING_CONSTANT * Vector3DotProduct(delta_velocity, delta_position) / current_length;
-
-    const Vector3 force_a = Vector3Scale(delta_position, -(hooke + dampening) / current_length);
-    const Vector3 force_b = Vector3Scale(force_a, -1.0);
-
-    _a.force = Vector3Add(_a.force, force_a);
-    _b.force = Vector3Add(_b.force, force_b);
+    positions[m] += delta_pos + acc;
+    previous_positions[m] = pos;
 }

 auto mass_spring_system::clear() -> void
 {
-    masses.clear();
+    positions.clear();
+    previous_positions.clear();
+    velocities.clear();
+    forces.clear();
    springs.clear();
    tree.nodes.clear();
 }
@ -77,24 +74,27 @@ auto mass_spring_system::add_mass() -> void
    // };
    // position = Vector3Scale(Vector3Normalize(position), REST_LENGTH * 2.0);

-    masses.emplace_back(Vector3Zero());
+    positions.emplace_back(Vector3Zero());
+    previous_positions.emplace_back(Vector3Zero());
+    velocities.emplace_back(Vector3Zero());
+    forces.emplace_back(Vector3Zero());
 }

 auto mass_spring_system::add_spring(size_t a, size_t b) -> void
 {
    // Update masses to be located along a random walk when adding the springs
-    const mass& mass_a = masses.at(a);
-    mass& mass_b = masses.at(b);
+    const Vector3& mass_a = positions[a];
+    const Vector3& mass_b = positions[b];

-    Vector3 offset{
-        static_cast<float>(GetRandomValue(-100, 100)), static_cast<float>(GetRandomValue(-100, 100)),
-        static_cast<float>(GetRandomValue(-100, 100))
-    };
+    Vector3 offset{static_cast<float>(GetRandomValue(-100, 100)),
+                   static_cast<float>(GetRandomValue(-100, 100)),
+                   static_cast<float>(GetRandomValue(-100, 100))};
    offset = Vector3Normalize(offset) * REST_LENGTH;

    // If the offset moves the mass closer to the current center of mass, flip it
    if (!tree.nodes.empty()) {
-        const Vector3 mass_center_direction = Vector3Subtract(mass_a.position, tree.nodes.at(0).mass_center);
+        const Vector3 mass_center_direction =
+            Vector3Subtract(positions[a], tree.nodes[0].mass_center);
        const float mass_center_distance = Vector3Length(mass_center_direction);

        if (mass_center_distance > 0 && Vector3DotProduct(offset, mass_center_direction) < 0.0f) {
@ -102,8 +102,8 @@ auto mass_spring_system::add_spring(size_t a, size_t b) -> void
        }
    }

-    mass_b.position = mass_a.position + offset;
-    mass_b.previous_position = mass_b.position;
+    positions[b] = mass_a + offset;
+    previous_positions[b] = mass_b;

    // infoln("Adding spring: ({}, {}, {})->({}, {}, {})", mass_a.position.x, mass_a.position.y,
    // mass_a.position.z,
@ -114,118 +114,92 @@ auto mass_spring_system::add_spring(size_t a, size_t b) -> void

 auto mass_spring_system::clear_forces() -> void
 {
-    #ifdef TRACY
+#ifdef TRACY
    ZoneScoped;
-    #endif
-
-    for (auto& m : masses) {
-        m.clear_force();
-    }
-}
-
-auto mass_spring_system::calculate_spring_forces() -> void
-{
-    #ifdef TRACY
-    ZoneScoped;
-    #endif
-
-    for (const auto s : springs) {
-        mass& a = masses.at(s.a);
-        mass& b = masses.at(s.b);
-        spring::calculate_spring_force(a, b);
-    }
-}
-
-#ifdef THREADPOOL
-auto mass_spring_system::set_thread_name(size_t idx) -> void
-{
-    BS::this_thread::set_os_thread_name(std::format("bh-worker-{}", idx));
-}
 #endif

-auto mass_spring_system::build_octree() -> void
+    memset(forces.data(), 0, forces.size() * sizeof(Vector3));
+}
+
+auto mass_spring_system::calculate_spring_forces(
+    const std::optional<BS::thread_pool<>* const> thread_pool) -> void
 {
-    #ifdef TRACY
+#ifdef TRACY
    ZoneScoped;
-    #endif
+#endif

-    tree.nodes.clear();
-    tree.nodes.reserve(masses.size() * 2);
+    const auto solve_spring_force = [&](const int i) { calculate_spring_force(i); };

-    // Compute bounding box around all masses
-    Vector3 min{FLT_MAX, FLT_MAX, FLT_MAX};
-    Vector3 max{-FLT_MAX, -FLT_MAX, -FLT_MAX};
-    for (const auto& m : masses) {
-        min.x = std::min(min.x, m.position.x);
-        max.x = std::max(max.x, m.position.x);
-        min.y = std::min(min.y, m.position.y);
-        max.y = std::max(max.y, m.position.y);
-        min.z = std::min(min.z, m.position.z);
-        max.z = std::max(max.z, m.position.z);
-    }
-
-    // Pad the bounding box
-    constexpr float pad = 1.0;
-    min = Vector3Subtract(min, Vector3Scale(Vector3One(), pad));
-    max = Vector3Add(max, Vector3Scale(Vector3One(), pad));
-
-    // Make it cubic (so subdivisions are balanced)
-    const float max_extent = std::max({max.x - min.x, max.y - min.y, max.z - min.z});
-    max = Vector3Add(min, Vector3Scale(Vector3One(), max_extent));
-
-    // Root node spans the entire area
-    const int root = tree.create_empty_leaf(min, max);
-
-    for (size_t i = 0; i < masses.size(); ++i) {
-        tree.insert(root, static_cast<int>(i), masses[i].position, MASS, 0);
+    if (thread_pool) {
+        (*thread_pool)
+            ->submit_loop(0, springs.size(), solve_spring_force, SMALL_TASK_BLOCK_SIZE)
+            .wait();
+    } else {
+        for (size_t i = 0; i < springs.size(); ++i) {
+            solve_spring_force(i);
+        }
    }
 }

-auto mass_spring_system::calculate_repulsion_forces() -> void
+auto mass_spring_system::calculate_repulsion_forces(
+    const std::optional<BS::thread_pool<>* const> thread_pool) -> void
 {
-    #ifdef TRACY
+#ifdef TRACY
    ZoneScoped;
-    #endif
+#endif

-    build_octree();
-
-    auto solve_octree = [&](const int i)
+    const auto solve_octree = [&](const int i)
    {
-        const Vector3 force = tree.calculate_force(0, masses[i].position);
-        masses[i].force = Vector3Add(masses[i].force, force);
+        const Vector3 force = tree.calculate_force(0, positions[i]);
+        forces[i] += force;
    };

    // Calculate forces using Barnes-Hut
-    #ifdef THREADPOOL
-    const BS::multi_future<void> loop_future = threads.submit_loop(0, masses.size(), solve_octree, 256);
-    loop_future.wait();
-    #else
-    for (size_t i = 0; i < masses.size(); ++i) {
-        solve_octree(i);
+    if (thread_pool) {
+        (*thread_pool)
+            ->submit_loop(0, positions.size(), solve_octree, LARGE_TASK_BLOCK_SIZE)
+            .wait();
+    } else {
+        for (size_t i = 0; i < positions.size(); ++i) {
+            solve_octree(i);
+        }
    }
-    #endif
 }

-auto mass_spring_system::verlet_update(const float delta_time) -> void
+auto mass_spring_system::update(const float dt,
+                                const std::optional<BS::thread_pool<>* const> thread_pool) -> void
 {
-    #ifdef TRACY
+#ifdef TRACY
    ZoneScoped;
-    #endif
+#endif

-    for (auto& m : masses) {
-        m.verlet_update(delta_time);
+    const auto update = [&](const int i) { verlet_update(i, dt); };
+
+    if (thread_pool) {
+        (*thread_pool)->submit_loop(0, positions.size(), update, SMALL_TASK_BLOCK_SIZE).wait();
+    } else {
+        for (size_t i = 0; i < positions.size(); ++i) {
+            update(i);
+        }
    }
 }

-auto mass_spring_system::center_masses() -> void
+auto mass_spring_system::center_masses(const std::optional<BS::thread_pool<>* const> thread_pool)
+    -> void
 {
    Vector3 mean = Vector3Zero();
-    for (const auto& m : masses) {
-        mean += m.position;
+    for (const Vector3& pos : positions) {
+        mean += pos;
    }
-    mean /= static_cast<float>(masses.size());
+    mean /= static_cast<float>(positions.size());

-    for (auto& m : masses) {
-        m.position -= mean;
+    const auto center_mass = [&](const int i) { positions[i] -= mean; };
+
+    if (thread_pool) {
+        (*thread_pool)->submit_loop(0, positions.size(), center_mass, SMALL_TASK_BLOCK_SIZE).wait();
+    } else {
+        for (size_t i = 0; i < positions.size(); ++i) {
+            center_mass(i);
+        }
    }
-}
+}
--- a/src/octree.cpp
+++ b/src/octree.cpp
@ -1,13 +1,9 @@
 #include "octree.hpp"
 #include "config.hpp"
-#include "util.hpp"

+#include <cfloat>
 #include <raymath.h>

-#ifdef TRACY
-#include <tracy/Tracy.hpp>
-#endif
-
 auto octree::node::child_count() const -> int
 {
    int child_count = 0;
@ -19,16 +15,6 @@ auto octree::node::child_count() const -> int
    return child_count;
 }

-auto octree::create_empty_leaf(const Vector3& box_min, const Vector3& box_max) -> int
-{
-    node n;
-    n.box_min = box_min;
-    n.box_max = box_max;
-    nodes.emplace_back(n);
-
-    return static_cast<int>(nodes.size() - 1);
-}
-
 auto octree::get_octant(const int node_idx, const Vector3& pos) const -> int
 {
    const node& n = nodes[node_idx];
@ -75,6 +61,16 @@ auto octree::get_child_bounds(const int node_idx, const int octant) const -> std
    return std::make_pair(min, max);
 }

+auto octree::create_empty_leaf(const Vector3& box_min, const Vector3& box_max) -> int
+{
+    node n;
+    n.box_min = box_min;
+    n.box_max = box_max;
+    nodes.emplace_back(n);
+
+    return static_cast<int>(nodes.size() - 1);
+}
+
 auto octree::insert(const int node_idx, const int mass_id, const Vector3& pos, const float mass,
                    const int depth) -> void
 {
@ -153,6 +149,44 @@ auto octree::insert(const int node_idx, const int mass_id, const Vector3& pos, c
    nodes[node_idx].mass_total = new_mass;
 }

+auto octree::build_octree(octree& t, const std::vector<Vector3>& positions) -> void
+{
+    #ifdef TRACY
+    ZoneScoped;
+    #endif
+
+    t.nodes.clear();
+    t.nodes.reserve(positions.size() * 2);
+
+    // Compute bounding box around all masses
+    Vector3 min{FLT_MAX, FLT_MAX, FLT_MAX};
+    Vector3 max{-FLT_MAX, -FLT_MAX, -FLT_MAX};
+    for (const auto& [x, y, z] : positions) {
+        min.x = std::min(min.x, x);
+        max.x = std::max(max.x, x);
+        min.y = std::min(min.y, y);
+        max.y = std::max(max.y, y);
+        min.z = std::min(min.z, z);
+        max.z = std::max(max.z, z);
+    }
+
+    // Pad the bounding box
+    constexpr float pad = 1.0;
+    min = Vector3Subtract(min, Vector3Scale(Vector3One(), pad));
+    max = Vector3Add(max, Vector3Scale(Vector3One(), pad));
+
+    // Make it cubic (so subdivisions are balanced)
+    const float max_extent = std::max({max.x - min.x, max.y - min.y, max.z - min.z});
+    max = Vector3Add(min, Vector3Scale(Vector3One(), max_extent));
+
+    // Root node spans the entire area
+    const int root = t.create_empty_leaf(min, max);
+
+    for (size_t i = 0; i < positions.size(); ++i) {
+        t.insert(root, static_cast<int>(i), positions[i], MASS, 0);
+    }
+}
+
 auto octree::calculate_force(const int node_idx, const Vector3& pos) const -> Vector3
 {
    if (node_idx < 0) {
--- a/src/orbit_camera.cpp
+++ b/src/orbit_camera.cpp
@ -4,10 +4,6 @@
 #include <raylib.h>
 #include <raymath.h>

-#ifdef TRACY
-    #include <tracy/Tracy.hpp>
-#endif
-
 auto orbit_camera::rotate(const Vector2 last_mouse, const Vector2 mouse) -> void
 {
    const auto [dx, dy] = Vector2Subtract(mouse, last_mouse);
@ -73,4 +69,4 @@ auto orbit_camera::update(const Vector3& current_target, const Vector3& mass_cen
    camera.target = target;
    camera.fovy = fov;
    camera.projection = projection;
-}
+}
--- a/src/puzzle.cpp
+++ b/src/puzzle.cpp
@ -3,70 +3,16 @@
 #include <algorithm>
 #include <boost/unordered/unordered_flat_map.hpp>

-auto puzzle::block::create_repr(const uint8_t x, const uint8_t y, const uint8_t w, const uint8_t h, const bool t,
+auto puzzle::block::create_repr(const uint8_t x,
+                                const uint8_t y,
+                                const uint8_t w,
+                                const uint8_t h,
+                                const bool t,
                                const bool i) -> uint16_t
 {
    return block().set_x(x).set_y(y).set_width(w).set_height(h).set_target(t).set_immovable(i).repr & ~INVALID;
 }

-auto puzzle::block::set_x(const uint8_t x) const -> block
-{
-    if (x > 7) {
-        throw std::invalid_argument("Block x-position out of bounds");
-    }
-
-    block b = *this;
-    set_bits(b.repr, X_S, X_E, x);
-    return b;
-}
-
-auto puzzle::block::set_y(const uint8_t y) const -> block
-{
-    if (y > 7) {
-        throw std::invalid_argument("Block y-position out of bounds");
-    }
-
-    block b = *this;
-    set_bits(b.repr, Y_S, Y_E, y);
-    return b;
-}
-
-auto puzzle::block::set_width(const uint8_t width) const -> block
-{
-    if (width - 1 > 7) {
-        throw std::invalid_argument("Block width out of bounds");
-    }
-
-    block b = *this;
-    set_bits(b.repr, WIDTH_S, WIDTH_E, width - 1u);
-    return b;
-}
-
-auto puzzle::block::set_height(const uint8_t height) const -> block
-{
-    if (height - 1 > 7) {
-        throw std::invalid_argument("Block height out of bounds");
-    }
-
-    block b = *this;
-    set_bits(b.repr, HEIGHT_S, HEIGHT_E, height - 1u);
-    return b;
-}
-
-auto puzzle::block::set_target(const bool target) const -> block
-{
-    block b = *this;
-    set_bits(b.repr, TARGET_S, TARGET_E, target);
-    return b;
-}
-
-auto puzzle::block::set_immovable(const bool immovable) const -> block
-{
-    block b = *this;
-    set_bits(b.repr, IMMOVABLE_S, IMMOVABLE_E, immovable);
-    return b;
-}
-
 auto puzzle::block::unpack_repr() const -> std::tuple<uint8_t, uint8_t, uint8_t, uint8_t, bool, bool>
 {
    const uint8_t x = get_x();
@ -79,46 +25,29 @@ auto puzzle::block::unpack_repr() const -> std::tuple<uint8_t, uint8_t, uint8_t,
    return {x, y, w, h, t, i};
 }

-auto puzzle::block::get_x() const -> uint8_t
+auto puzzle::block::hash() const -> size_t
 {
-    return get_bits(repr, X_S, X_E);
+    return std::hash<uint16_t>{}(repr);
 }

-auto puzzle::block::get_y() const -> uint8_t
+auto puzzle::block::position_independent_hash() const -> size_t
 {
-    return get_bits(repr, Y_S, Y_E);
-}
-
-auto puzzle::block::get_width() const -> uint8_t
-{
-    return get_bits(repr, WIDTH_S, WIDTH_E) + 1u;
-}
-
-auto puzzle::block::get_height() const -> uint8_t
-{
-    return get_bits(repr, HEIGHT_S, HEIGHT_E) + 1u;
-}
-
-auto puzzle::block::get_target() const -> bool
-{
-    return get_bits(repr, TARGET_S, TARGET_E);
-}
-
-auto puzzle::block::get_immovable() const -> bool
-{
-    return get_bits(repr, IMMOVABLE_S, IMMOVABLE_E);
+    uint16_t r = repr;
+    clear_bits(r, X_S, X_E);
+    clear_bits(r, Y_S, Y_E);
+    return std::hash<uint16_t>{}(r);
 }

 auto puzzle::block::valid() const -> bool
 {
-    const auto [x, y, w, h, t, i] = unpack_repr();
-
-    if (t && i) {
+    // This means the first bit is set, marking the block as empty
+    if (repr & INVALID) {
        return false;
    }

-    // This means the first bit is set, marking the block as empty
-    if (repr & INVALID) {
+    const auto [x, y, w, h, t, i] = unpack_repr();
+
+    if (t && i) {
        return false;
    }

@ -171,8 +100,13 @@ auto puzzle::create_meta(const std::tuple<uint8_t, uint8_t, uint8_t, uint8_t, bo
    return m;
 }

-auto puzzle::create_repr(const uint8_t w, const uint8_t h, const uint8_t tx, const uint8_t ty, const bool r,
-                         const bool g, const std::array<uint16_t, MAX_BLOCKS>& b) -> repr_cooked
+auto puzzle::create_repr(const uint8_t w,
+                         const uint8_t h,
+                         const uint8_t tx,
+                         const uint8_t ty,
+                         const bool r,
+                         const bool g,
+                         const std::array<uint16_t, MAX_BLOCKS>& b) -> repr_cooked
 {
    repr_cooked repr = puzzle().set_width(w).set_height(h).set_goal_x(tx).set_goal_y(ty).set_restricted(r).set_goal(g).
                                set_blocks(b).repr.cooked;
@ -180,7 +114,9 @@ auto puzzle::create_repr(const uint8_t w, const uint8_t h, const uint8_t tx, con
    return repr;
 }

-auto puzzle::create_repr(const uint64_t byte_0, const uint64_t byte_1, const uint64_t byte_2,
+auto puzzle::create_repr(const uint64_t byte_0,
+                         const uint64_t byte_1,
+                         const uint64_t byte_2,
                         const uint64_t byte_3) -> repr_cooked
 {
    repr_u repr{};
@ -199,64 +135,6 @@ auto puzzle::create_repr(const std::string& string_repr) -> repr_cooked
    return *repr;
 }

-auto puzzle::set_restricted(const bool restricted) const -> puzzle
-{
-    uint16_t meta = repr.cooked.meta;
-    set_bits(meta, RESTRICTED_S, RESTRICTED_E, restricted);
-    return puzzle(meta, repr.cooked.blocks);
-}
-
-auto puzzle::set_width(const uint8_t width) const -> puzzle
-{
-    if (width - 1 > MAX_WIDTH) {
-        throw "Board width out of bounds";
-    }
-
-    uint16_t meta = repr.cooked.meta;
-    set_bits(meta, WIDTH_S, WIDTH_E, width - 1u);
-    return puzzle(meta, repr.cooked.blocks);
-}
-
-auto puzzle::set_height(const uint8_t height) const -> puzzle
-{
-    if (height - 1 > MAX_HEIGHT) {
-        throw "Board height out of bounds";
-    }
-
-    uint16_t meta = repr.cooked.meta;
-    set_bits(meta, HEIGHT_S, HEIGHT_E, height - 1u);
-    return puzzle(meta, repr.cooked.blocks);
-}
-
-auto puzzle::set_goal(const bool goal) const -> puzzle
-{
-    uint16_t meta = repr.cooked.meta;
-    set_bits(meta, GOAL_S, GOAL_E, goal);
-    return puzzle(meta, repr.cooked.blocks);
-}
-
-auto puzzle::set_goal_x(const uint8_t target_x) const -> puzzle
-{
-    if (target_x >= MAX_WIDTH) {
-        throw "Board target x out of bounds";
-    }
-
-    uint16_t meta = repr.cooked.meta;
-    set_bits(meta, GOAL_X_S, GOAL_X_E, target_x);
-    return puzzle(meta, repr.cooked.blocks);
-}
-
-auto puzzle::set_goal_y(const uint8_t target_y) const -> puzzle
-{
-    if (target_y >= MAX_HEIGHT) {
-        throw "Board target y out of bounds";
-    }
-
-    uint16_t meta = repr.cooked.meta;
-    set_bits(meta, GOAL_Y_S, GOAL_Y_E, target_y);
-    return puzzle(meta, repr.cooked.blocks);
-}
-
 auto puzzle::set_blocks(std::array<uint16_t, MAX_BLOCKS> blocks) const -> puzzle
 {
    puzzle p = *this;
@ -277,36 +155,6 @@ auto puzzle::unpack_meta() const -> std::tuple<uint8_t, uint8_t, uint8_t, uint8_
    return {w, h, tx, ty, r, g};
 }

-auto puzzle::get_restricted() const -> bool
-{
-    return get_bits(repr.cooked.meta, RESTRICTED_S, RESTRICTED_E);
-}
-
-auto puzzle::get_width() const -> uint8_t
-{
-    return get_bits(repr.cooked.meta, WIDTH_S, WIDTH_E) + 1u;
-}
-
-auto puzzle::get_height() const -> uint8_t
-{
-    return get_bits(repr.cooked.meta, HEIGHT_S, HEIGHT_E) + 1u;
-}
-
-auto puzzle::get_goal() const -> bool
-{
-    return get_bits(repr.cooked.meta, GOAL_S, GOAL_E);
-}
-
-auto puzzle::get_goal_x() const -> uint8_t
-{
-    return get_bits(repr.cooked.meta, GOAL_X_S, GOAL_X_E);
-}
-
-auto puzzle::get_goal_y() const -> uint8_t
-{
-    return get_bits(repr.cooked.meta, GOAL_Y_S, GOAL_Y_E);
-}
-
 auto puzzle::hash() const -> size_t
 {
    size_t h = 0;
@ -542,7 +390,7 @@ auto puzzle::try_get_invalid_reason() const -> std::optional<std::string>

    const auto [w, h, gx, gy, r, g] = unpack_meta();

-    infoln("Validating puzzle {}", string_repr());
+    // traceln("Validating puzzle \"{}\"", string_repr());

    const std::optional<block>& b = try_get_target_block();
    if (get_goal() && !b) {
@ -896,64 +744,8 @@ auto puzzle::try_move_block_at(const uint8_t x, const uint8_t y, const direction
    return p;
 }

-auto puzzle::try_move_block_at_fast(const uint64_t bitmap, const uint8_t block_idx,
-                                    const direction dir) const -> std::optional<puzzle>
-{
-    const block b = block(repr.cooked.blocks[block_idx]);
-    const auto [bx, by, bw, bh, bt, bi] = b.unpack_repr();
-    if (bi) {
-        return std::nullopt;
-    }
-
-    const auto [w, h, gx, gy, r, g] = unpack_meta();
-    const int dirs = r ? b.principal_dirs() : nor | eas | sou | wes;
-
-    // Get target block
-    int _target_x = bx;
-    int _target_y = by;
-    switch (dir) {
-    case nor:
-        if (!(dirs & nor) || _target_y < 1) {
-            return std::nullopt;
-        }
-        --_target_y;
-        break;
-    case eas:
-        if (!(dirs & eas) || _target_x + bw >= w) {
-            return std::nullopt;
-        }
-        ++_target_x;
-        break;
-    case sou:
-        if (!(dirs & sou) || _target_y + bh >= h) {
-            return std::nullopt;
-        }
-        ++_target_y;
-        break;
-    case wes:
-        if (!(dirs & wes) || _target_x < 1) {
-            return std::nullopt;
-        }
-        --_target_x;
-        break;
-    }
-
-    // Check collisions
-    const uint64_t bm = bitmap_clear_block(bitmap, b);
-    if (bitmap_check_collision(bm, b, dir)) {
-        return std::nullopt;
-    }
-
-    // Replace block
-    const std::array<uint16_t, MAX_BLOCKS> blocks = sorted_replace(repr.cooked.blocks, block_idx,
-                                                                   block::create_repr(
-                                                                       _target_x, _target_y, bw, bh, bt));
-
-    // This constructor doesn't sort
-    return puzzle(std::make_tuple(w, h, gx, gy, r, g), blocks);
-}
-
-auto puzzle::sorted_replace(std::array<uint16_t, MAX_BLOCKS> blocks, const uint8_t idx,
+auto puzzle::sorted_replace(std::array<uint16_t, MAX_BLOCKS> blocks,
+                            const uint8_t idx,
                            const uint16_t new_val) -> std::array<uint16_t, MAX_BLOCKS>
 {
    // Remove old entry
@ -987,95 +779,69 @@ auto puzzle::blocks_bitmap() const -> uint64_t
        }

        auto [x, y, w, h, t, im] = b.unpack_repr();
+        const uint8_t width = get_width();

        for (int dy = 0; dy < h; ++dy) {
            for (int dx = 0; dx < w; ++dx) {
-                bitmap |= 1ULL << ((y + dy) * 8 + (x + dx));
+                bitmap_set_bit(bitmap, width, x + dx, y + dy);
            }
        }
    }
    return bitmap;
 }

-inline auto puzzle::bitmap_set_bit(const uint64_t bitmap, const uint8_t x, const uint8_t y) -> uint64_t
+auto puzzle::blocks_bitmap_h() const -> uint64_t
 {
-    return bitmap & ~(1ULL << (y * 8 + x));
-}
+    uint64_t bitmap = 0;
+    for (uint8_t i = 0; i < MAX_BLOCKS; ++i) {
+        block b(repr.cooked.blocks[i]);
+        if (!b.valid()) {
+            break;
+        }
+        const int dirs = b.principal_dirs();
+        if (!(dirs & eas)) {
+            continue;
+        }

-inline auto puzzle::bitmap_get_bit(const uint64_t bitmap, const uint8_t x, const uint8_t y) -> bool
-{
-    return bitmap & (1ULL << (y * 8 + x));
-}
+        auto [x, y, w, h, t, im] = b.unpack_repr();
+        const uint8_t width = get_width();

-auto puzzle::bitmap_clear_block(uint64_t bitmap, const block b) -> uint64_t
-{
-    const auto [x, y, w, h, t, i] = b.unpack_repr();
-
-    for (int dy = 0; dy < h; ++dy) {
-        for (int dx = 0; dx < w; ++dx) {
-            bitmap = bitmap_set_bit(bitmap, x + dx, y + dy);
+        for (int dy = 0; dy < h; ++dy) {
+            for (int dx = 0; dx < w; ++dx) {
+                bitmap_set_bit(bitmap, width, x + dx, y + dy);
+            }
        }
    }
-
    return bitmap;
 }

-auto puzzle::bitmap_check_collision(const uint64_t bitmap, const block b) -> bool
+auto puzzle::blocks_bitmap_v() const -> uint64_t
 {
-    const auto [x, y, w, h, t, i] = b.unpack_repr();
+    uint64_t bitmap = 0;
+    for (uint8_t i = 0; i < MAX_BLOCKS; ++i) {
+        block b(repr.cooked.blocks[i]);
+        if (!b.valid()) {
+            break;
+        }
+        const int dirs = b.principal_dirs();
+        if (!(dirs & sou)) {
+            continue;
+        }

-    for (int dy = 0; dy < h; ++dy) {
-        for (int dx = 0; dx < w; ++dx) {
-            if (bitmap_get_bit(bitmap, x + dx, y + dy)) {
-                return true; // collision
+        auto [x, y, w, h, t, im] = b.unpack_repr();
+        const uint8_t width = get_width();
+
+        for (int dy = 0; dy < h; ++dy) {
+            for (int dx = 0; dx < w; ++dx) {
+                bitmap_set_bit(bitmap, width, x + dx, y + dy);
            }
        }
    }
-
-    return false;
-}
-
-auto puzzle::bitmap_check_collision(const uint64_t bitmap, const block b, const direction dir) -> bool
-{
-    const auto [x, y, w, h, t, i] = b.unpack_repr();
-
-    switch (dir) {
-    case nor: // Check the row above: (x...x+w-1, y-1)
-        for (int dx = 0; dx < w; ++dx) {
-            if (bitmap_get_bit(bitmap, x + dx, y - 1)) {
-                return true;
-            }
-        }
-        break;
-    case sou: // Check the row below: (x...x+w-1, y+h)
-        for (int dx = 0; dx < w; ++dx) {
-            if (bitmap_get_bit(bitmap, x + dx, y + h)) {
-                return true;
-            }
-        }
-        break;
-    case wes: // Check the column left: (x-1, y...y+h-1)
-        for (int dy = 0; dy < h; ++dy) {
-            if (bitmap_get_bit(bitmap, x - 1, y + dy)) {
-                return true;
-            }
-        }
-        break;
-    case eas: // Check the column right: (x+w, y...y+h-1)
-        for (int dy = 0; dy < h; ++dy) {
-            if (bitmap_get_bit(bitmap, x + w, y + dy)) {
-                return true;
-            }
-        }
-        break;
-    }
-    return false;
+    return bitmap;
 }

 auto puzzle::explore_state_space() const -> std::pair<std::vector<puzzle>, std::vector<std::pair<size_t, size_t>>>
 {
-    const std::chrono::high_resolution_clock::time_point start = std::chrono::high_resolution_clock::now();
-
    std::vector<puzzle> state_pool;
    boost::unordered_flat_map<puzzle, std::size_t, puzzle_hasher> state_indices;
    std::vector<std::pair<size_t, size_t>> links;
@ -1083,6 +849,11 @@ auto puzzle::explore_state_space() const -> std::pair<std::vector<puzzle>, std::
    // Buffer for all states we want to call GetNextStates() on
    std::vector<size_t> queue; // indices into state_pool

+    #ifdef WIP
+    // Store an index to the blocks array of a state for each occupied bitmap cell
+    std::array<uint8_t, 64> bitmap_block_indices;
+    #endif
+
    // Start with the current state
    state_indices.emplace(*this, 0);
    state_pool.push_back(*this);
@ -1095,6 +866,24 @@ auto puzzle::explore_state_space() const -> std::pair<std::vector<puzzle>, std::
        // Make a copy because references might be invalidated when inserting into the vector
        const puzzle current = state_pool[current_idx];

+        #ifdef WIP
+        // Build bitmap-block indices
+        for (size_t i = 0; i < MAX_BLOCKS; ++i) {
+            const block b = block(current.repr.cooked.blocks[i]);
+            const auto [bx, by, bw, bh, bt, bi] = b.unpack_repr();
+            if (!b.valid()) {
+                break;
+            }
+
+            for (uint8_t x = bx; x < bx + bw; ++x) {
+                for (uint8_t y = by; y < by + bh; ++y) {
+                    bitmap_block_indices[y * current.get_width() + x] = i;
+                }
+            }
+        }
+        #endif
+
+        // TODO: I can just dispatch different functions depending on if the board is restricted or contains walls
        current.for_each_adjacent([&](const puzzle& p)
        {
            auto [it, inserted] = state_indices.emplace(p, state_pool.size());
@ -1106,10 +895,239 @@ auto puzzle::explore_state_space() const -> std::pair<std::vector<puzzle>, std::
        });
    }

-    const std::chrono::high_resolution_clock::time_point end = std::chrono::high_resolution_clock::now();
-
-    infoln("Explored puzzle. Took {} ms.", std::chrono::duration_cast<std::chrono::milliseconds>(end - start).count());
-    infoln("State space has size {} with {} transitions.", state_pool.size(), links.size());
-
    return {std::move(state_pool), std::move(links)};
+}
+
+auto puzzle::get_cluster_id_and_solution() const -> std::pair<puzzle, bool>
+{
+    const auto& [puzzles, moves] = explore_state_space();
+    bool solution = false;
+    puzzle min = puzzles[0];
+    for (size_t i = 0; i < puzzles.size(); ++i) {
+        if (puzzles[i] < min) {
+            min = puzzles[i];
+        }
+        if (puzzles[i].goal_reached()) {
+            solution = true;
+        }
+    }
+    return {min, solution};
+}
+
+auto puzzle::bitmap_find_first_empty(const uint64_t bitmap, int& x, int& y) const -> bool
+{
+    x = 0;
+    y = 0;
+
+    // Bitmap is empty of first slot is empty
+    if (bitmap_is_empty(bitmap) || !(bitmap & 1u)) {
+        return true;
+    }
+
+    // Bitmap is full
+    if (bitmap_is_full(bitmap)) {
+        return false;
+    }
+
+    // Find the next more significant empty bit (we know the first slot is full)
+    int ls_set = 0;
+    bool next_set = true;
+    while (next_set && ls_set < get_width() * get_height() - 1) {
+        next_set = bitmap & (1ul << (ls_set + 1));
+        ++ls_set;
+    }
+
+    x = ls_set % get_width();
+    y = ls_set / get_width();
+    return true;
+}
+
+auto puzzle::generate_block_sequences(
+    const boost::unordered_flat_set<block, block_hasher2, block_equal2>& permitted_blocks,
+    const block target_block,
+    const size_t max_blocks,
+    std::vector<block>& current_sequence,
+    const int current_area,
+    const int board_area,
+    const std::function<void(const std::vector<block>&)>& callback) -> void
+{
+    if (!current_sequence.empty()) {
+        callback(current_sequence);
+    }
+
+    if (current_sequence.size() == max_blocks) {
+        return;
+    }
+
+    for (const block b : permitted_blocks) {
+        const int new_area = current_area + b.get_width() * b.get_height();
+        if (new_area > board_area) {
+            continue;
+        }
+
+        // Explore all sequences with the block placed, then continue the loop
+        current_sequence.push_back(b);
+        generate_block_sequences(permitted_blocks,
+                                 target_block,
+                                 max_blocks,
+                                 current_sequence,
+                                 new_area,
+                                 board_area,
+                                 callback);
+        current_sequence.pop_back();
+    }
+}
+
+auto puzzle::place_block_sequence(const puzzle& p,
+                                  const uint64_t& bitmap,
+                                  const std::tuple<uint8_t, uint8_t, uint8_t, uint8_t, bool, bool>& p_repr,
+                                  const std::vector<block>& sequence,
+                                  const block target_block,
+                                  const std::tuple<uint8_t, uint8_t, uint8_t, uint8_t>& target_block_pos_range,
+                                  const bool has_target,
+                                  const size_t index,
+                                  const std::function<void(const puzzle&)>& callback) -> void
+{
+    if (index == sequence.size()) {
+        // All blocks placed
+        callback(p);
+        return;
+    }
+
+    if (!has_target && p.get_restricted()) {
+        // Place target block (restricted movement)
+        const auto [txs, tys, txe, tye] = target_block_pos_range;
+        for (int tx = txs; tx <= txe; ++tx) {
+            for (int ty = tys; ty <= tye; ++ty) {
+                block t = target_block;
+                t = t.set_x(tx);
+                t = t.set_y(ty);
+
+                if (!p.covers(t)) {
+                    continue;
+                }
+
+                const std::array<uint16_t, MAX_BLOCKS> blocks = sorted_replace(p.repr.cooked.blocks, 0, t.repr);
+                const puzzle next_p = puzzle(p_repr, blocks);
+
+                uint64_t next_bm = bitmap;
+                next_p.bitmap_set_block(next_bm, t);
+
+                // Place the remaining blocks for each possible target block configuration
+                // traceln("Generating block sequence for target at {},{}", tx, ty);
+                next_p.place_block_sequence(next_p,
+                                            next_bm,
+                                            p_repr,
+                                            sequence,
+                                            target_block,
+                                            target_block_pos_range,
+                                            true,
+                                            index,
+                                            callback);
+            }
+        }
+
+        return;
+    }
+
+
+    if (!has_target && !p.get_restricted()) {
+        // Place target block (free movement)
+
+        // TODO
+    }
+
+    int x, y;
+    if (!p.bitmap_find_first_empty(bitmap, x, y)) {
+        // No space remaining
+        callback(p);
+        return;
+    }
+
+    block b = sequence[index];
+    b = b.set_x(static_cast<uint8_t>(x));
+    b = b.set_y(static_cast<uint8_t>(y));
+
+    // Place the next block and call the resulting subtree, then remove the block and continue here
+    if (!p.bitmap_check_collision(bitmap, b) && p.covers(b)) {
+        // Shift the sequence by 1 (index + 1), because the target block is inserted separately
+        const std::array<uint16_t, MAX_BLOCKS> blocks = sorted_replace(p.repr.cooked.blocks, index + 1, b.repr);
+        const puzzle next_p = puzzle(p_repr, blocks);
+
+        uint64_t next_bm = bitmap;
+        next_p.bitmap_set_block(next_bm, b);
+
+        next_p.place_block_sequence(next_p,
+                                    next_bm,
+                                    p_repr,
+                                    sequence,
+                                    target_block,
+                                    target_block_pos_range,
+                                    true,
+                                    index + 1,
+                                    callback);
+    }
+
+    // Create an empty cell and call the resulting subtree (without advancing the block index)
+    uint64_t next_bm = bitmap;
+    bitmap_set_bit(next_bm, p.get_width(), b.get_x(), b.get_y());
+    p.place_block_sequence(p, next_bm, p_repr, sequence, target_block, target_block_pos_range, true, index, callback);
+}
+
+auto puzzle::explore_puzzle_space(const boost::unordered_flat_set<block, block_hasher2, block_equal2>& permitted_blocks,
+                                  const block target_block,
+                                  const std::tuple<uint8_t, uint8_t, uint8_t, uint8_t>& target_block_pos_range,
+                                  const size_t max_blocks,
+                                  const std::optional<BS::thread_pool<>* const> thread_pool) const ->
+    boost::unordered_flat_set<puzzle, puzzle_hasher>
+{
+    const auto [w, h, gx, gy, r, g] = unpack_meta();
+
+    // Implemented in the slowest, stupidest way for now:
+    // 1. Iterate through all possible permitted_blocks permutations using recursive tree descent
+    // 2. Find the cluster id of the permutation by populating the entire state space
+    //    - We could do some preprocessing to quickly reduce the numeric value
+    //      of the state and check if its already contained in visited_clusters,
+    //      this could save some state space calculations.
+    // 3. Add it to visited_clusters if unseen
+
+    std::mutex mtx;
+
+    boost::unordered_flat_set<puzzle, puzzle_hasher> visited_clusters;
+
+    // TODO: Can't even parallelize this. Or just start at different initial puzzles?
+    const puzzle empty_puzzle = puzzle(w, h, gx, gy, r, g);
+    const auto board_repr = std::make_tuple(w, h, gx, gy, r, g);
+    std::vector<block> current_sequence;
+    int total = 0;
+    generate_block_sequences(permitted_blocks,
+                             target_block,
+                             max_blocks - 1, // Make space for the target block
+                             current_sequence,
+                             target_block.get_width() * target_block.get_height(), // Starting area
+                             get_width() * get_height(),
+                             [&](const std::vector<block>& sequence)
+                             {
+                                 place_block_sequence(empty_puzzle,
+                                                      0,
+                                                      board_repr,
+                                                      sequence,
+                                                      target_block,
+                                                      target_block_pos_range,
+                                                      false,
+                                                      0,
+                                                      [&](const puzzle& p)
+                                                      {
+                                                          const auto [cluster_id, winnable] = p.get_cluster_id_and_solution();
+                                                          std::lock_guard<std::mutex> lock(mtx);
+                                                          ++total;
+                                                          if (winnable) {
+                                                              visited_clusters.emplace(cluster_id);
+                                                          }
+                                                      });
+                             });
+
+    infoln("Found {} of {} clusters with a solution", visited_clusters.size(), total);
+
+    return visited_clusters;
 }
--- a/src/renderer.cpp
+++ b/src/renderer.cpp
@ -5,33 +5,29 @@
 #include <raymath.h>
 #include <rlgl.h>

-#ifdef TRACY
-    #include <tracy/Tracy.hpp>
-#endif
-
 auto renderer::update_texture_sizes() -> void
 {
    if (!IsWindowResized()) {
        return;
    }

-    UnloadRenderTexture(render_target);
+    UnloadRenderTexture(graph_target);
    UnloadRenderTexture(klotski_target);
    UnloadRenderTexture(menu_target);

    const int width = GetScreenWidth() / 2;
    const int height = GetScreenHeight() - MENU_HEIGHT;

-    render_target = LoadRenderTexture(width, height);
+    graph_target = LoadRenderTexture(width, height);
    klotski_target = LoadRenderTexture(width, height);
    menu_target = LoadRenderTexture(width * 2, MENU_HEIGHT);
 }

 auto renderer::draw_mass_springs(const std::vector<Vector3>& masses) -> void
 {
-#ifdef TRACY
+    #ifdef TRACY
    ZoneScoped;
-#endif
+    #endif

    if (masses.size() != state.get_state_count()) {
        // Because the physics run in a different thread, it might need time to catch up
@ -40,16 +36,16 @@ auto renderer::draw_mass_springs(const std::vector<Vector3>& masses) -> void

    // Prepare connection batching
    {
-#ifdef TRACY
-        ZoneNamedN(prepare_masses, "PrepareConnectionsBatching", true);
-#endif
+        #ifdef TRACY
+        ZoneNamedN(prepare_connections, "PrepareConnectionsBatching", true);
+        #endif

        connections.clear();
        connections.reserve(state.get_target_count());
        if (input.connect_solutions) {
            for (const size_t& _state : state.get_winning_indices()) {
-                const Vector3& current_mass = masses.at(state.get_current_index());
-                const Vector3& winning_mass = masses.at(_state);
+                const Vector3& current_mass = masses[state.get_current_index()];
+                const Vector3& winning_mass = masses[_state];
                connections.emplace_back(current_mass, winning_mass);
                DrawLine3D(current_mass, winning_mass, Fade(TARGET_BLOCK_COLOR, 0.5));
            }
@ -58,9 +54,9 @@ auto renderer::draw_mass_springs(const std::vector<Vector3>& masses) -> void

    // Prepare cube instancing
    {
-#ifdef TRACY
+        #ifdef TRACY
        ZoneNamedN(prepare_masses, "PrepareMassInstancing", true);
-#endif
+        #endif

        if (masses.size() < DRAW_VERTICES_LIMIT) {
            // Don't have to reserve, capacity is already set to DRAW_VERTICES_LIMIT in constructor
@ -73,12 +69,10 @@ auto renderer::draw_mass_springs(const std::vector<Vector3>& masses) -> void

                // Normal vertex
                Color c = VERTEX_COLOR;
-                if ((input.mark_solutions || input.mark_path) &&
-                    state.get_winning_indices().contains(mass)) {
+                if ((input.mark_solutions || input.mark_path) && state.get_winning_indices().contains(mass)) {
                    // Winning vertex
                    c = VERTEX_TARGET_COLOR;
-                } else if ((input.mark_solutions || input.mark_path) &&
-                           state.get_path_indices().contains(mass)) {
+                } else if ((input.mark_solutions || input.mark_path) && state.get_path_indices().contains(mass)) {
                    // Path vertex
                    c = VERTEX_PATH_COLOR;
                } else if (mass == state.get_starting_index()) {
@ -98,21 +92,21 @@ auto renderer::draw_mass_springs(const std::vector<Vector3>& masses) -> void
        rlUpdateVertexBuffer(color_vbo_id, colors.data(), colors.size() * sizeof(Color), 0);
    }

-    BeginTextureMode(render_target);
+    BeginTextureMode(graph_target);
    ClearBackground(RAYWHITE);
    BeginMode3D(camera.camera);

    // Draw springs (batched)
    {
-#ifdef TRACY
+        #ifdef TRACY
        ZoneNamedN(draw_springs, "DrawSprings", true);
-#endif
+        #endif

        rlBegin(RL_LINES);
        for (const auto& [from, to] : state.get_links()) {
            if (masses.size() > from && masses.size() > to) {
-                const auto& [ax, ay, az] = masses.at(from);
-                const auto& [bx, by, bz] = masses.at(to);
+                const auto& [ax, ay, az] = masses[from];
+                const auto& [bx, by, bz] = masses[to];
                rlColor4ub(EDGE_COLOR.r, EDGE_COLOR.g, EDGE_COLOR.b, EDGE_COLOR.a);
                rlVertex3f(ax, ay, az);
                rlVertex3f(bx, by, bz);
@ -123,17 +117,16 @@ auto renderer::draw_mass_springs(const std::vector<Vector3>& masses) -> void

    // Draw masses (instanced)
    {
-#ifdef TRACY
+        #ifdef TRACY
        ZoneNamedN(draw_masses, "DrawMasses", true);
-#endif
+        #endif

        if (masses.size() < DRAW_VERTICES_LIMIT) {
            // NOTE: I don't know if drawing all this inside a shader would make it
            //       much faster... The amount of data sent to the GPU would be
            //       reduced (just positions instead of matrices), but is this
            //       noticable for < 100000 cubes?
-            DrawMeshInstanced(cube_instance, vertex_mat, transforms.data(),
-                              masses.size()); // NOLINT(*-narrowing-conversions)
+            DrawMeshInstanced(cube_instance, vertex_mat, transforms.data(), masses.size()); // NOLINT(*-narrowing-conversions)
        }
    }

@ -152,9 +145,8 @@ auto renderer::draw_mass_springs(const std::vector<Vector3>& masses) -> void
    // Mark current state
    const size_t current_index = state.get_current_index();
    if (masses.size() > current_index) {
-        const Vector3& current_mass = masses.at(current_index);
-        DrawCube(current_mass, VERTEX_SIZE * 2, VERTEX_SIZE * 2, VERTEX_SIZE * 2,
-                 VERTEX_CURRENT_COLOR);
+        const Vector3& current_mass = masses[current_index];
+        DrawCube(current_mass, VERTEX_SIZE * 2, VERTEX_SIZE * 2, VERTEX_SIZE * 2, VERTEX_CURRENT_COLOR);
    }

    EndMode3D();
@ -163,9 +155,9 @@ auto renderer::draw_mass_springs(const std::vector<Vector3>& masses) -> void

 auto renderer::draw_klotski() const -> void
 {
-#ifdef TRACY
+    #ifdef TRACY
    ZoneScoped;
-#endif
+    #endif

    BeginTextureMode(klotski_target);
    ClearBackground(RAYWHITE);
@ -177,9 +169,9 @@ auto renderer::draw_klotski() const -> void

 auto renderer::draw_menu() const -> void
 {
-#ifdef TRACY
+    #ifdef TRACY
    ZoneScoped;
-#endif
+    #endif

    BeginTextureMode(menu_target);
    ClearBackground(RAYWHITE);
@ -194,32 +186,28 @@ auto renderer::draw_textures(const int fps, const int ups, const size_t mass_cou
 {
    BeginDrawing();

-    DrawTextureRec(menu_target.texture,
-                   Rectangle(0, 0, menu_target.texture.width, -menu_target.texture.height),
+    DrawTextureRec(menu_target.texture, Rectangle(0, 0, menu_target.texture.width, -menu_target.texture.height),
                   Vector2(0, 0), WHITE);
    DrawTextureRec(klotski_target.texture,
                   Rectangle(0, 0, klotski_target.texture.width, -klotski_target.texture.height),
                   Vector2(0, MENU_HEIGHT), WHITE);
-    DrawTextureRec(render_target.texture,
-                   Rectangle(0, 0, render_target.texture.width, -render_target.texture.height),
+    DrawTextureRec(graph_target.texture, Rectangle(0, 0, graph_target.texture.width, -graph_target.texture.height),
                   Vector2(GetScreenWidth() / 2.0f, MENU_HEIGHT), WHITE);

    // Draw borders
    DrawRectangleLinesEx(Rectangle(0, 0, GetScreenWidth(), MENU_HEIGHT), 1.0f, BLACK);
-    DrawRectangleLinesEx(
-        Rectangle(0, MENU_HEIGHT, GetScreenWidth() / 2.0f, GetScreenHeight() - MENU_HEIGHT), 1.0f,
-        BLACK);
+    DrawRectangleLinesEx(Rectangle(0, MENU_HEIGHT, GetScreenWidth() / 2.0f, GetScreenHeight() - MENU_HEIGHT), 1.0f,
+                         BLACK);
    DrawRectangleLinesEx(Rectangle(GetScreenWidth() / 2.0f, MENU_HEIGHT, GetScreenWidth() / 2.0f,
-                                   GetScreenHeight() - MENU_HEIGHT),
-                         1.0f, BLACK);
+                                   GetScreenHeight() - MENU_HEIGHT), 1.0f, BLACK);

    gui.draw(fps, ups, mass_count, spring_count);

    EndDrawing();
 }

-auto renderer::render(const std::vector<Vector3>& masses, const int fps, const int ups,
-                      const size_t mass_count, const size_t spring_count) -> void
+auto renderer::render(const std::vector<Vector3>& masses, const int fps, const int ups, const size_t mass_count,
+                      const size_t spring_count) -> void
 {
    update_texture_sizes();

@ -227,4 +215,4 @@ auto renderer::render(const std::vector<Vector3>& masses, const int fps, const i
    draw_klotski();
    draw_menu();
    draw_textures(fps, ups, mass_count, spring_count);
-}
+}
--- a/src/state_manager.cpp
+++ b/src/state_manager.cpp
@ -2,17 +2,10 @@
 #include "graph_distances.hpp"
 #include "util.hpp"

-#include <fstream>
-#include <ios>
-
-#ifdef TRACY
-#include <tracy/Tracy.hpp>
-#endif
-
 auto state_manager::synced_try_insert_state(const puzzle& state) -> size_t
 {
    if (state_indices.contains(state)) {
-        return state_indices.at(state);
+        return state_indices[state];
    }

    const size_t index = state_pool.size();
@ -77,78 +70,27 @@ auto state_manager::synced_clear_statespace() -> void
    physics.clear_cmd();
 }

-auto state_manager::parse_preset_file(const std::string& _preset_file) -> bool
+auto state_manager::save_current_to_preset_file(const std::string& preset_comment) -> void
 {
-    preset_file = _preset_file;
-
-    std::ifstream file(preset_file);
-    if (!file) {
-        infoln("Preset file \"{}\" couldn't be loaded.", preset_file);
-        return false;
+    if (append_preset_file(preset_file, preset_comment, get_current_state())) {
+        current_preset = preset_states.size();
+        reload_preset_file();
    }
-
-    std::string line;
-    std::vector<std::string> comment_lines;
-    std::vector<std::string> preset_lines;
-    while (std::getline(file, line)) {
-        if (line.starts_with("S")) {
-            preset_lines.push_back(line);
-        } else if (line.starts_with("#")) {
-            comment_lines.push_back(line);
-        }
-    }
-
-    if (preset_lines.empty() || comment_lines.size() != preset_lines.size()) {
-        infoln("Preset file \"{}\" couldn't be loaded.", preset_file);
-        return false;
-    }
-
-    preset_states.clear();
-    for (const auto& preset : preset_lines) {
-        // Each char is a bit
-        const puzzle& p = puzzle(preset);
-
-        if (const std::optional<std::string>& reason = p.try_get_invalid_reason()) {
-            preset_states = {puzzle(4, 5, 0, 0, true, false)};
-            infoln("Preset file \"{}\" contained invalid presets: {}", preset_file, *reason);
-            return false;
-        }
-        preset_states.emplace_back(p);
-    }
-    preset_comments = comment_lines;
-
-    infoln("Loaded {} presets from \"{}\".", preset_lines.size(), preset_file);
-
-    return true;
 }

-auto state_manager::append_preset_file(const std::string& preset_name) -> bool
+auto state_manager::reload_preset_file() -> void
 {
-    infoln(R"(Saving preset "{}" to "{}")", preset_name, preset_file);
-
-    if (get_current_state().try_get_invalid_reason()) {
-        return false;
+    const auto [presets, comments] = parse_preset_file(preset_file);
+    if (!presets.empty()) {
+        preset_states = presets;
+        preset_comments = comments;
    }
-
-    std::ofstream file(preset_file, std::ios_base::app | std::ios_base::out);
-    if (!file) {
-        infoln("Preset file \"{}\" couldn't be loaded.", preset_file);
-        return false;
-    }
-
-    file << "\n# " << preset_name << "\n" << get_current_state().string_repr() << std::flush;
-
-    infoln("Refreshing presets...");
-    if (parse_preset_file(preset_file)) {
-        load_preset(preset_states.size() - 1);
-    }
-
-    return true;
+    load_preset(current_preset);
 }

 auto state_manager::load_preset(const size_t preset) -> void
 {
-    clear_graph_and_add_current(preset_states.at(preset));
+    clear_graph_and_add_current(preset_states[preset]);
    current_preset = preset;
    edited = false;
 }
@ -270,8 +212,8 @@ auto state_manager::goto_most_distant_state() -> void
    int max_distance = 0;
    size_t max_distance_index = 0;
    for (size_t i = 0; i < node_target_distances.distances.size(); ++i) {
-        if (node_target_distances.distances.at(i) > max_distance) {
-            max_distance = node_target_distances.distances.at(i);
+        if (node_target_distances.distances[i] > max_distance) {
+            max_distance = node_target_distances.distances[i];
            max_distance_index = i;
        }
    }
@ -285,7 +227,7 @@ auto state_manager::goto_closest_target_state() -> void
        return;
    }

-    update_current_state(get_state(node_target_distances.nearest_targets.at(current_state_index)));
+    update_current_state(get_state(node_target_distances.nearest_targets[current_state_index]));
 }

 auto state_manager::populate_graph() -> void
@ -298,17 +240,23 @@ auto state_manager::populate_graph() -> void
    const puzzle s = get_starting_state();
    const puzzle p = get_current_state();

-
    // Clear the graph first so we don't add duplicates somehow
    synced_clear_statespace();

    // Explore the entire statespace starting from the current state
+    const std::chrono::high_resolution_clock::time_point start = std::chrono::high_resolution_clock::now();
+
    const auto& [states, _links] = s.explore_state_space();
    synced_insert_statespace(states, _links);

-    current_state_index = state_indices.at(p);
+    const std::chrono::high_resolution_clock::time_point end = std::chrono::high_resolution_clock::now();
+
+    infoln("Explored puzzle. Took {}ms.", std::chrono::duration_cast<std::chrono::milliseconds>(end - start).count());
+    infoln("State space has size {} with {} transitions.", state_pool.size(), links.size());
+
+    current_state_index = state_indices[p];
    previous_state_index = current_state_index;
-    starting_state_index = state_indices.at(s);
+    starting_state_index = state_indices[s];

    // Search for cool stuff
    populate_winning_indices();
@ -393,7 +341,7 @@ auto state_manager::get_starting_index() const -> size_t

 auto state_manager::get_state(const size_t index) const -> const puzzle&
 {
-    return state_pool.at(index);
+    return state_pool[index];
 }

 auto state_manager::get_current_state() const -> const puzzle&
@ -468,7 +416,7 @@ auto state_manager::get_preset_count() const -> size_t

 auto state_manager::get_current_preset_comment() const -> const std::string&
 {
-    return preset_comments.at(current_preset);
+    return preset_comments[current_preset];
 }

 auto state_manager::has_history() const -> bool
--- a/src/threaded_physics.cpp
+++ b/src/threaded_physics.cpp
@ -1,6 +1,7 @@
 #include "threaded_physics.hpp"
 #include "config.hpp"
 #include "mass_spring_system.hpp"
+#include "util.hpp"

 #include <chrono>
 #include <raylib.h>
@ -8,20 +9,30 @@
 #include <utility>
 #include <vector>

-#ifdef TRACY
-#include <tracy/Tracy.hpp>
+#ifdef ASYNC_OCTREE
+auto threaded_physics::set_octree_pool_thread_name(size_t idx) -> void
+{
+    BS::this_thread::set_os_thread_name(std::format("octree-{}", idx));
+    // traceln("Using thread \"{}\"", BS::this_thread::get_os_thread_name().value_or("INVALID NAME"));
+}
 #endif

-auto threaded_physics::physics_thread(physics_state& state) -> void
+auto threaded_physics::physics_thread(physics_state& state, const std::optional<BS::thread_pool<>* const> thread_pool) -> void
 {
-    #ifdef THREADPOOL
-    BS::this_thread::set_os_thread_name("physics");
-    #endif
-
    mass_spring_system mass_springs;

+    #ifdef ASYNC_OCTREE
+    BS::this_thread::set_os_thread_name("physics");
+
+    BS::thread_pool<> octree_thread(1, set_octree_pool_thread_name);
+    std::future<void> octree_future;
+    octree tree_buffer;
+    size_t last_mass_count = 0;
+    infoln("Using asynchronous octree builder.");
+    #endif
+
    const auto visitor = overloads{
-        [&](const struct add_mass& am)
+        [&](const struct add_mass&)
        {
            mass_springs.add_mass();
        },
@ -29,7 +40,7 @@ auto threaded_physics::physics_thread(physics_state& state) -> void
        {
            mass_springs.add_spring(as.a, as.b);
        },
-        [&](const struct clear_graph& cg)
+        [&](const struct clear_graph&)
        {
            mass_springs.clear();
        },
@ -66,26 +77,58 @@ auto threaded_physics::physics_thread(physics_state& state) -> void
            }
        }

-        if (mass_springs.masses.empty()) {
+        if (mass_springs.positions.empty()) {
            std::this_thread::sleep_for(std::chrono::milliseconds(1));
            continue;
        }

        // Physics update
        if (physics_accumulator.count() > TIMESTEP) {
+            #ifdef ASYNC_OCTREE
+            // Snapshot the positions so mass_springs is not mutating the vector while the octree is building
+            std::vector<Vector3> positions = mass_springs.positions;
+
+            // Start building the octree for the next physics update.
+            // Move the snapshot into the closure so it doesn't get captured by reference (don't use [&])
+            octree_future = octree_thread.submit_task([&tree_buffer, positions = std::move(positions)]()
+            {
+                octree::build_octree(tree_buffer, positions);
+            });
+
+            // Rebuild the tree synchronously if we changed the number of masses to not use
+            // an empty tree from the last frame in the frame where the graph was generated
+            if (last_mass_count != mass_springs.positions.size()) {
+                traceln("Rebuilding octree synchronously because graph size changed");
+                octree::build_octree(mass_springs.tree, mass_springs.positions);
+                last_mass_count = mass_springs.positions.size();
+            }
+            #else
+            octree::build_octree(mass_springs.tree, mass_springs.positions);
+            #endif
+
            mass_springs.clear_forces();
-            mass_springs.calculate_spring_forces();
-            mass_springs.calculate_repulsion_forces();
-            mass_springs.verlet_update(TIMESTEP * SIM_SPEED);
+            mass_springs.calculate_spring_forces(thread_pool);
+            mass_springs.calculate_repulsion_forces(thread_pool);
+            mass_springs.update(TIMESTEP * SIM_SPEED, thread_pool);

            // This is only helpful if we're drawing a grid at (0, 0, 0). Otherwise, it's just
            // expensive and yields no benefit since we can lock the camera to the center of mass
-            // cheaply. mass_springs.center_masses();
+            // cheaply.
+            // mass_springs.center_masses(thread_pool);

            ++loop_iterations;
            physics_accumulator -= std::chrono::duration<double>(TIMESTEP);
        }

+        #ifdef ASYNC_OCTREE
+        // Wait for the octree to be built
+        if (octree_future.valid()) {
+            octree_future.wait();
+            octree_future = std::future<void>{};
+            std::swap(mass_springs.tree, tree_buffer);
+        }
+        #endif
+
        // Publish the positions for the renderer (copy)
        #ifdef TRACY
        FrameMarkStart("PhysicsThreadProduceLock");
@ -114,16 +157,16 @@ auto threaded_physics::physics_thread(physics_state& state) -> void
            if (mass_springs.tree.nodes.empty()) {
                state.mass_center = Vector3Zero();
            } else {
-                state.mass_center = mass_springs.tree.nodes.at(0).mass_center;
+                state.mass_center = mass_springs.tree.nodes[0].mass_center;
            }

            state.masses.clear();
-            state.masses.reserve(mass_springs.masses.size());
-            for (const auto& mass : mass_springs.masses) {
-                state.masses.emplace_back(mass.position);
+            state.masses.reserve(mass_springs.positions.size());
+            for (const Vector3& pos : mass_springs.positions) {
+                state.masses.emplace_back(pos);
            }

-            state.mass_count = mass_springs.masses.size();
+            state.mass_count = mass_springs.positions.size();
            state.spring_count = mass_springs.springs.size();

            state.data_ready = true;
@ -131,14 +174,25 @@ auto threaded_physics::physics_thread(physics_state& state) -> void
        }
        // Notify the rendering thread that new data is available
        state.data_ready_cnd.notify_all();
-        #ifdef TRACY
-        FrameMarkEnd("PhysicsThreadProduceLock");

-        FrameMarkEnd("PhysicsThread");
+        #ifdef TRACY
+        FrameMarkEnd("PhysicsThreadProduceLock"); FrameMarkEnd("PhysicsThread");
        #endif
    }
 }

+auto threaded_physics::clear_cmd() -> void
+{
+    {
+        #ifdef TRACY
+        std::lock_guard<LockableBase(std::mutex)> lock(state.command_mtx);
+        #else
+        std::lock_guard<std::mutex> lock(state.command_mtx);
+        #endif
+        state.pending_commands.emplace(clear_graph{});
+    }
+}
+
 auto threaded_physics::add_mass_cmd() -> void
 {
    {
@ -163,18 +217,6 @@ auto threaded_physics::add_spring_cmd(const size_t a, const size_t b) -> void
    }
 }

-auto threaded_physics::clear_cmd() -> void
-{
-    {
-        #ifdef TRACY
-        std::lock_guard<LockableBase(std::mutex)> lock(state.command_mtx);
-        #else
-        std::lock_guard<std::mutex> lock(state.command_mtx);
-        #endif
-        state.pending_commands.emplace(clear_graph{});
-    }
-}
-
 auto threaded_physics::add_mass_springs_cmd(const size_t num_masses,
                                            const std::vector<std::pair<size_t, size_t>>& springs) -> void
 {
--- a/src/user_interface.cpp
+++ b/src/user_interface.cpp
@ -7,10 +7,6 @@
 #define RAYGUI_IMPLEMENTATION
 #include <raygui.h>

-#ifdef TRACY
-#include <tracy/Tracy.hpp>
-#endif
-
 auto user_interface::grid::update_bounds(const int _x, const int _y, const int _width, const int _height,
                                         const int _columns, const int _rows) -> void
 {
@ -639,13 +635,13 @@ auto user_interface::draw_save_preset_popup() -> void

    // Returns the pressed button index
    const int button = GuiTextInputBox(popup_bounds(), "Save as Preset", "Enter Preset Name", "Ok;Cancel",
-                                       preset_name.data(), 255, nullptr);
+                                       preset_comment.data(), 255, nullptr);
    if (button == 1) {
-        state.append_preset_file(preset_name.data());
+        state.save_current_to_preset_file(preset_comment.data());
    }
    if (button == 0 || button == 1 || button == 2) {
        save_window = false;
-        TextCopy(preset_name.data(), "\0");
+        TextCopy(preset_comment.data(), "\0");
    }
 }

--- a/test/bitmap.cpp
+++ b/test/bitmap.cpp
@ -0,0 +1,74 @@
+// ReSharper disable CppLocalVariableMayBeConst
+#include "puzzle.hpp"
+
+#include <random>
+#include <catch2/catch_test_macros.hpp>
+
+TEST_CASE("bitmap_is_full all bits set", "[puzzle][board]")
+{
+    puzzle p1(5, 5);
+    puzzle p2(3, 4);
+    puzzle p3(5, 4);
+    puzzle p4(3, 7);
+    uint64_t bitmap = -1;
+
+    REQUIRE(p1.bitmap_is_full(bitmap));
+    REQUIRE(p2.bitmap_is_full(bitmap));
+    REQUIRE(p3.bitmap_is_full(bitmap));
+    REQUIRE(p4.bitmap_is_full(bitmap));
+}
+
+TEST_CASE("bitmap_is_full no bits set", "[puzzle][board]")
+{
+    puzzle p1(5, 5);
+    puzzle p2(3, 4);
+    puzzle p3(5, 4);
+    puzzle p4(3, 7);
+    uint64_t bitmap = 0;
+
+    REQUIRE_FALSE(p1.bitmap_is_full(bitmap));
+    REQUIRE_FALSE(p2.bitmap_is_full(bitmap));
+    REQUIRE_FALSE(p3.bitmap_is_full(bitmap));
+    REQUIRE_FALSE(p4.bitmap_is_full(bitmap));
+}
+
+TEST_CASE("bitmap_is_full necessary bits set", "[puzzle][board]")
+{
+    puzzle p1(5, 5);
+    puzzle p2(3, 4);
+    puzzle p3(5, 4);
+    puzzle p4(3, 7);
+
+    uint64_t bitmap1 = (1ull << 25) - 1; // 5 * 5
+    uint64_t bitmap2 = (1ull << 12) - 1;  // 3 * 4
+    uint64_t bitmap3 = (1ull << 20) - 1; // 5 * 4
+    uint64_t bitmap4 = (1ull << 21) - 1; // 3 * 7
+
+    REQUIRE(p1.bitmap_is_full(bitmap1));
+    REQUIRE(p2.bitmap_is_full(bitmap2));
+    REQUIRE(p3.bitmap_is_full(bitmap3));
+    REQUIRE(p4.bitmap_is_full(bitmap4));
+}
+
+TEST_CASE("bitmap_is_full necessary bits not set", "[puzzle][board]")
+{
+    puzzle p1(5, 5);
+    puzzle p2(3, 4);
+    puzzle p3(5, 4);
+    puzzle p4(3, 7);
+
+    uint64_t bitmap1 = (1ull << 25) - 1; // 5 * 5
+    uint64_t bitmap2 = (1ull << 12) - 1;  // 3 * 4
+    uint64_t bitmap3 = (1ull << 20) - 1; // 5 * 4
+    uint64_t bitmap4 = (1ull << 21) - 1; // 3 * 7
+
+    bitmap1 &= ~(1ull << 12);
+    bitmap2 &= ~(1ull << 6);
+    bitmap3 &= ~(1ull << 8);
+    bitmap4 &= ~(1ull << 18);
+
+    REQUIRE_FALSE(p1.bitmap_is_full(bitmap1));
+    REQUIRE_FALSE(p2.bitmap_is_full(bitmap2));
+    REQUIRE_FALSE(p3.bitmap_is_full(bitmap3));
+    REQUIRE_FALSE(p4.bitmap_is_full(bitmap4));
+}
--- a/test/bitmap_find_first_empty.cpp
+++ b/test/bitmap_find_first_empty.cpp
@ -0,0 +1,266 @@
+// ReSharper disable CppLocalVariableMayBeConst
+#include "puzzle.hpp"
+
+#include <random>
+#include <catch2/catch_test_macros.hpp>
+#include <catch2/generators/catch_generators.hpp>
+
+static auto board_mask(const int w, const int h) -> uint64_t
+{
+    const int cells = w * h;
+    if (cells == 64) {
+        return ~0ULL;
+    }
+    return (1ULL << cells) - 1ULL;
+}
+
+TEST_CASE("Empty board returns (0,0)", "[puzzle][board]")
+{
+    puzzle p(5, 5);
+
+    int x = -1, y = -1;
+    REQUIRE(p.bitmap_find_first_empty(0ULL, x, y));
+
+    REQUIRE(x == 0);
+    REQUIRE(y == 0);
+}
+
+TEST_CASE("Full board detection respects width*height only", "[puzzle][board]")
+{
+    auto [w, h] = GENERATE(std::tuple{3, 3}, std::tuple{4, 4}, std::tuple{5, 6}, std::tuple{8, 8});
+
+    puzzle p(w, h);
+
+    uint64_t mask = board_mask(w, h);
+
+    int x = -1, y = -1;
+
+    REQUIRE_FALSE(p.bitmap_find_first_empty(mask, x, y));
+
+    // Bits outside board should not affect fullness
+    REQUIRE_FALSE(p.bitmap_find_first_empty(mask | (~mask), x, y));
+}
+
+TEST_CASE("Single empty cell at various positions", "[puzzle][board]")
+{
+    auto [w, h] = GENERATE(std::tuple{3, 3}, std::tuple{4, 4}, std::tuple{5, 5}, std::tuple{8, 8});
+
+    puzzle p(w, h);
+
+    int cells = w * h;
+
+    auto empty_index = GENERATE_COPY(values<int>({ 0, cells / 2, cells - 1}));
+
+    uint64_t bitmap = board_mask(w, h);
+    bitmap &= ~(1ULL << empty_index);
+
+    int x = -1, y = -1;
+    REQUIRE(p.bitmap_find_first_empty(bitmap, x, y));
+
+    REQUIRE(x == empty_index % w);
+    REQUIRE(y == empty_index / w);
+}
+
+TEST_CASE("Bits outside board are ignored", "[puzzle][board]")
+{
+    puzzle p(3, 3); // 9 cells
+
+    uint64_t mask = board_mask(3, 3);
+
+    // Board is full
+    uint64_t bitmap = mask;
+
+    // Set extra bits outside 9 cells
+    bitmap |= (1ULL << 20);
+    bitmap |= (1ULL << 63);
+
+    int x = -1, y = -1;
+    REQUIRE_FALSE(p.bitmap_find_first_empty(bitmap, x, y));
+}
+
+TEST_CASE("First empty found in forward search branch", "[puzzle][branch]")
+{
+    puzzle p(4, 4); // 16 cells
+
+    // Only MSB (within board) set
+    uint64_t bitmap = (1ULL << 15);
+
+    int x = -1, y = -1;
+    REQUIRE(p.bitmap_find_first_empty(bitmap, x, y));
+
+    REQUIRE(x == 0);
+    REQUIRE(y == 0);
+}
+
+TEST_CASE("Backward search branch finds gap before MSB cluster", "[puzzle][branch]")
+{
+    puzzle p(4, 4); // 16 cells
+
+    // Set bits 15,14,13 but leave 12 empty
+    uint64_t bitmap = (1ULL << 15) | (1ULL << 14) | (1ULL << 13);
+
+    int x = -1, y = -1;
+    REQUIRE(p.bitmap_find_first_empty(bitmap, x, y));
+
+    REQUIRE(x == 0);
+    REQUIRE(y == 0);
+}
+
+TEST_CASE("Rectangular board coordinate mapping", "[puzzle][rect]")
+{
+    puzzle p(5, 3); // 15 cells
+
+    int empty_index = 11;
+
+    uint64_t bitmap = board_mask(5, 3);
+    bitmap &= ~(1ULL << empty_index);
+
+    int x = -1, y = -1;
+    REQUIRE(p.bitmap_find_first_empty(bitmap, x, y));
+
+    REQUIRE(x == empty_index % 5);
+    REQUIRE(y == empty_index / 5);
+}
+
+TEST_CASE("Non-64-sized board near limit", "[puzzle][limit]")
+{
+    puzzle p(7, 8); // 56 cells
+
+    uint64_t mask = board_mask(7, 8);
+
+    // Full board should return false
+    int x = -1, y = -1;
+    REQUIRE_FALSE(p.bitmap_find_first_empty(mask, x, y));
+
+    // Clear highest valid cell
+    int empty_index = 55;
+    mask &= ~(1ULL << empty_index);
+
+    REQUIRE(p.bitmap_find_first_empty(mask, x, y));
+    REQUIRE(x == empty_index % 7);
+    REQUIRE(y == empty_index / 7);
+}
+
+// --- Oracle: find first zero bit inside board ---
+static auto oracle_find_first_empty(uint64_t bitmap, int w, int h, int& x, int& y) -> bool
+{
+    int cells = w * h;
+
+    for (int i = 0; i < cells; ++i) {
+        if ((bitmap & (1ULL << i)) == 0) {
+            x = i % w;
+            y = i / w;
+            return true;
+        }
+    }
+
+    return false;
+}
+
+TEST_CASE("Oracle validation across board sizes 3x3 to 8x8", "[puzzle][oracle]")
+{
+    auto [w, h] = GENERATE(std::tuple{3, 3}, std::tuple{4, 4}, std::tuple{5, 5}, std::tuple{6, 6}, std::tuple{7, 7},
+                           std::tuple{8, 8}, std::tuple{3, 5}, std::tuple{5, 3}, std::tuple{7, 8}, std::tuple{8, 7});
+
+    puzzle p(w, h);
+
+    uint64_t mask = board_mask(w, h);
+
+    std::mt19937_64 rng(12345);
+    std::uniform_int_distribution<uint64_t> dist(0, UINT64_MAX);
+
+    for (int iteration = 0; iteration < 200; ++iteration) {
+        uint64_t bitmap = dist(rng);
+
+        int ox = -1, oy = -1;
+        bool oracle_result = oracle_find_first_empty(bitmap, w, h, ox, oy);
+
+        int x = -1, y = -1;
+        bool result = p.bitmap_find_first_empty(bitmap, x, y);
+
+        REQUIRE(result == oracle_result);
+
+        if (result) {
+            REQUIRE(x == ox);
+            REQUIRE(y == oy);
+        }
+    }
+}
+
+TEST_CASE("Bits set outside board only behaves as empty board", "[puzzle][outside]")
+{
+    puzzle p(3, 3); // 9 cells
+
+    uint64_t bitmap = (1ULL << 40) | (1ULL << 63);
+
+    int x = -1, y = -1;
+    REQUIRE(p.bitmap_find_first_empty(bitmap, x, y));
+
+    REQUIRE(x == 0);
+    REQUIRE(y == 0);
+}
+
+TEST_CASE("Last valid cell empty stresses upper bound", "[puzzle][boundary]")
+{
+    auto [w, h] = GENERATE(std::tuple{4, 4}, std::tuple{5, 6}, std::tuple{7, 8}, std::tuple{8, 8});
+
+    puzzle p(w, h);
+
+    int cells = w * h;
+    uint64_t bitmap = board_mask(w, h);
+
+    // Clear last valid bit
+    bitmap &= ~(1ULL << (cells - 1));
+
+    int x = -1, y = -1;
+    REQUIRE(p.bitmap_find_first_empty(bitmap, x, y));
+
+    REQUIRE(x == (cells - 1) % w);
+    REQUIRE(y == (cells - 1) / w);
+}
+
+TEST_CASE("Board sizes between 33 and 63 cells", "[puzzle][midrange]")
+{
+    auto [w, h] = GENERATE(std::tuple{6, 6}, // 36
+                           std::tuple{7, 6}, // 42
+                           std::tuple{7, 7}, // 49
+                           std::tuple{8, 7}, // 56
+                           std::tuple{7, 8}  // 56
+    );
+
+    puzzle p(w, h);
+
+    int cells = w * h;
+
+    for (int index : {31, 32, cells - 2}) {
+        if (index >= cells) continue;
+
+        uint64_t bitmap = board_mask(w, h);
+        bitmap &= ~(1ULL << index);
+
+        int x = -1, y = -1;
+        REQUIRE(p.bitmap_find_first_empty(bitmap, x, y));
+
+        REQUIRE(x == index % w);
+        REQUIRE(y == index / w);
+    }
+}
+
+TEST_CASE("Multiple holes choose lowest index", "[puzzle][multiple]")
+{
+    puzzle p(5, 5);
+
+    uint64_t bitmap = board_mask(5, 5);
+
+    // Clear several positions
+    bitmap &= ~(1ULL << 3);
+    bitmap &= ~(1ULL << 7);
+    bitmap &= ~(1ULL << 12);
+
+    int x = -1, y = -1;
+    REQUIRE(p.bitmap_find_first_empty(bitmap, x, y));
+
+    // Oracle expectation: index 3
+    REQUIRE(x == 3 % 5);
+    REQUIRE(y == 3 / 5);
+}
--- a/test/puzzle.cpp
+++ b/test/puzzle.cpp
Author	SHA1	Message	Date
Christoph Urlacher	d4f83e11db	add .desktop icon to package	2026-03-05 02:05:02 +01:00
Christoph Urlacher	db588bd57b	fix conditional threadpool include	2026-03-05 02:04:48 +01:00
Christoph Urlacher	49e5ed6906	fix windows build	2026-03-04 21:39:39 +01:00
Christoph Urlacher	08352dd997	don't pass a reference to a temporary to physics_thread	2026-03-04 21:31:01 +01:00
Christoph Urlacher	4e5ca6be6c	fix nix build error	2026-03-04 20:47:34 +01:00
Christoph Urlacher	a9c102298a	enable/disable threadpool from cmake	2026-03-04 20:45:25 +01:00
Christoph Urlacher	cc2aee3af4	dummy commit	2026-03-04 20:33:29 +01:00
Christoph Urlacher	e0f128f693	dummy commit	2026-03-04 20:29:33 +01:00
Christoph Urlacher	3b6919944c	add some abbrs to flake + wrap package to set the working dir	2026-03-04 20:28:38 +01:00
Christoph Urlacher	c9915852db	implement very slow puzzle space exploration	2026-03-04 20:23:16 +01:00
Christoph Urlacher	2d111f58da	add single state space benchmark + some tests	2026-03-04 19:08:00 +01:00
Christoph Urlacher	0a2788c1b4	update flake	2026-03-04 19:07:40 +01:00
Christoph Urlacher	7a5013295e	fix cross compilation for windows (disable boost.stacktrace + use sqrt instead of rsqrt)	2026-03-02 20:09:49 +01:00
Christoph Urlacher	9c48954a78	build octree async and reuse tree from last frame (disabled as it breaks physics)	2026-03-02 18:52:06 +01:00
Christoph Urlacher	d62d5c78bf	fix bug where interactive movement was using the fast move implementation instead of the interactive one	2026-03-02 14:38:03 +01:00
Christoph Urlacher	2a5f1b2ffd	update default camera settings	2026-03-02 14:37:41 +01:00
Christoph Urlacher	2ef2a29601	squash merge efficient-puzzle into main	2026-03-02 14:37:34 +01:00