nissy-core

h48.md (23837B)

---

 # The H48 optimal solver
 
 This document contains information on the H48 Rubik's Cube optimal solver.
 This solver is occasionally improved and optimized, and this document
 is updated accordingly to reflect the current implementation.
 
 I highly encourage the reader to check out Jaap Scherphuis'
 [Computer Puzzling page](https://www.jaapsch.net/puzzles/compcube.htm)
 before reading this document.
 
 Other backround material that will be referenced throughout the text includes:
 
 * [Cube coordinates](https://sebastiano.tronto.net/speedcubing/coordinates)
   by Sebastiano Tronto
 * [Nxopt](https://github.com/rokicki/cube20src/blob/master/nxopt.md)
   by Tomas Rokicki
 * [The Mathematics behind Cube Explorer](https://kociemba.org/cube.htm)
   by Herber Kociemba
 
 Initially I intended to give minimum information here, citing the
 above resources when needed, but I ended up rewriting even some of the
 basic things in this document.
 
 ## Optimal solvers
 
 ### IDA*
 
 The basic idea behind the solver is an iterative-deepening
 depth-first search with pruning tables, also known as IDA*. You
 can find a detailed explanation of this idea in
 [Jaap's page](https://www.jaapsch.net/puzzles/compcube.htm#tree).
 
 To summarize, the IDA* algorithm finds the shortest solution(s) for a
 scrambled puzzle by using multiple depth-first searches at increasing
 depth. A breadth-first search is not applicable in this scenario because
 of the large amount of memory required.  At each depth, the puzzle's state
 is evaluated to obtain an estimated lower bound for how many moves are
 required to solve the puzzle; this is done by employing large pruning
 tables, as well as other techniques. If the lower bound exceeds the
 number of moves required to reach the current target depth, the branch is
 pruned. Otherwise, every possible move is tried, except those that would
 "cancel" with the preceding moves, obtaining new positions at depth N+1.
 
 ### Pruning tables and coordinates
 
 A pruning table associates to a cube position a value that is a
 lower bound for the number of moves it takes to solve that position.
 To acces this value, the cube must be turned into an index for the
 given table. This is done using a
 [*coordinate*](https://sebastiano.tronto.net/speedcubing/coordinates),
 by which we mean a function from the set of (admissible, solvable) cube
 positions to the set of non-negative integers. The maximum value N of a
 coordinate must be less than the size of the pruning table; preferably
 there are no "gaps", that is the coordinate reaches all values between
 0 and N and the pruning table has size N+1.
 
 Coordinates are often derived from the cosets of a *target subgroup*.
 As an example, consider the pruning table that has one entry for each
 possible corner position, and the stored values denote the length of
 an optimal corner solution for that position. In this situation, the
 target subgroup is the set of all cube positions with solved corners
 (for the Mathematician: this is indeed a subgroup of the group of all
 cube positions). This group consists of `12! * 2^10` positions.  The set
 of cosets of this subgroup can be identified with the set of corner
 configurations, and it has size `8! * 3^7`.  By looking at corners,
 from each cube position we can compute a number from `0` to `8! * 3^7 - 1`
 and use it as an index for the pruning table.
 
 The pruning table can be filled in many ways. A common way is starting
 with the coordinate of the solved cube (often 0) and visiting the set
 of all coordinates with a breadth-first search. To do this, one needs
 to apply a cube move to a coordinate; this can either be done directly
 (for example using transistion tables, as explained in Jaap's page) or
 by converting the coordinate back into a valid cube position, applying
 the move and then computing the coordinate of the new position.
 
 The largest a pruning table is, the more information it contains,
 the faster the solvers that uses it is. It is therefore convenient to
 compress the tables as much as possible. This can be achieved by using
 4 bits per position to store the lower bound, or as little as 2 bits
 per position with some caveats (see Jaap's page or Rokicki's nxopt
 description for two possible ways to achieve this).  Tables that use
 only 1 bit per position are possible, but the only information they
 can provide is wether the given position requires more or less than a
 certain number of moves to be solved.
 
 ### Symmetries
 
 Another, extremely effective way of reducing the size of a pruning table
 is using symmetries. In the example above of the corner table, one can
 notice that all corner positions that are one quarter-turn away from being
 solved are essentially the same: they can be turned one into the other
 by rotating the cube and optionally applying a mirroring transformation.
 Indeed, every position belongs to a set of up to 48 positions that are
 "essentially the same" in this way.  Keeping a separate entry in the
 pruning table for each of these positions is a waste of memory.
 
 There are ways to reduce a coordinate by symmetry, grouping every
 transformation-equivalent position into the same index. You can
 find a description in my
 [cube coordinates page](https://sebastiano.tronto.net/speedcubing/coordinates).
 By doing so, the size of the pruning table is reduced without loss of
 information.
 
 Some well-known solvers (Cube Explorer, nxopt) do not take advantage
 of the full group of 48 symmetries of the cube, but they use only the
 16 symmetries that fix the UD axis. However, they make up for this by
 doing 3 pruning table lookups, one for each axis of the cube.
 
 Some cube positions are self-symmetric: when applying certain
 transformations to them, they remain the same. For example, the
 cube obtained by applying U2 to a solved cube is invariant with
 respect to the mirroring along the RL-plane. This fact has a couple
 of consequences:
 
 * Reducing by a group of N symmetries does not reduce the size of
   the pruning table by a factor of N, but by a slightly smaller
   factor. If the size of the coordinate is large, this fact is
   harmless.
 * When combining symmetry-reduced coordinates with other coordinates,
   one has to be extra careful with handling self-symmetric coordinates.
   An informal explanation of this phenomenon can be found in
   the [docs/nasty-symmetries.md document](./nasty-symmetries.md).
 
 ## The H48 solver
 
 ### The target subgroup: H48 coordinates
 
 The H48 solver uses a target group that is invariant under all 48
 symmetries. This group is defined as follows:
 
 * Each corner is placed in the *corner tetrad* it belongs to, and it
   is oriented. Here by corner tetrad we mean one of the two sets of
   corners {UFR, UBL, DFL, DBR} and {UFL, UBR, DFR, DBL}. For a corner
   that is placed in its own tetrad, the orientation does not depend on
   the reference axis chosen, and an oriented corner can be defined
   has having the U or D sticker facing U or D.
 * Each edge is placed in the *slice* it belongs to. The three edge slices
   are M = {UF, UB, DB, DF}, S = {UR, UL, DL, DR} and E = {FR, FL, BL, BR}.
 
 Other options are available for corners: for example, we could
 impose that the corner permutation is even or that corners are in
 *half-turn reduction* state, that is they are in the group generated
 by {U2, D2, R2, L2, F2, B2}. We settled on the corner group described
 above because it is easy to calculate and it gives enough options
 for pruning table sizes, depending on how edge orientation is
 considered (see the next section).
 
 We call the coordinate obtained from the target subgroup described
 above the **h0** coordinate. This coordinate has `C(8,4) * 3^7 *
 C(12,8) * C(8,4) = 5.304.568.500` possible values.  If the corner
 part of the coordinate is reduced by symmetry, it consists of only
 3393 elements, giving a total of `3393 * C(12,8) * C(8,4) =
 117.567.450` values. This is not very large for a pruning table,
 as with 2 bits per entry it would occupy less than 30MB. However,
 it is possible to take edge orientation (partially) into account
 to produce larger tables.
 
 ### Edge orientation
 
 Like for corners, the orientation of an edge is well-defined independently
 of any axis when said edge is placed in the slice it belongs to. We may
 modify the target subgroup defined in the previous section by imposing
 that the edges are all oriented. This yields a new coordinate, which we
 call **h11**, whose full size with symmetry-reduced corners is `3393 *
 C(12,8) * C(8,4) * 2^11 = 240.778.137.600`. A pruning table based on
 this coordinate with 2 bits per entry would occupy around 60GB, which
 is a little too much for most personal computers.
 
 One can wonder if it is possible to use a coordinate that considers
 the orientation of only *some* of the edges, which we may call **h1**
 to **h10**. Such coordinates do exist, but they are not invariant under
 the full symmetry group: indeed an edge whose orientation we keep track
 of could be moved to any of the untracked edges' positions by one of
 the symmetries, making the whole coordinate ill-defined.
 
 It is however possible to compute the symmetry-reduced pruning tables
 for these coordinates. One way to construct them is by taking the **h11**
 pruning table and "forgetting" about some of the edge orientation values,
 collapsing 2 (or a power thereof) values to one by taking the minimum.
 It is also possible to compute these tables directly, as explained in
 the **Pruning table computation** section below.
 
 ### Coordinate computation for pruning value estimation
 
 In order to access the pruning value for a given cube position, we first
 need to compute its coordinate value. Let's use for simplicity the **h0**
 coordinate, but everything will be valid for any other coordinate of
 the type described in the previous section.
 
 As described in the symmetric-composite coordinate section of
 [my coordinates page](https://sebastiano.tronto.net/speedcubing/coordinates),
 we first need to compute the part of the coordinate where the symmetry
 is applied, that is the corner separation + orientation (also called
 *cocsep*). The value of the coordinate depends on the equivalence
 class of the corner configuration, and is memorized in a table called
 `cocsepdata`. To avoid lengthy computations, the index for a cube position
 in this table is determined by the corner orientation (between `0` and
 `3^7-1`) and a binary representation of the corner separation part (an
 8-bit number with 4 zero digits and 4 one digits, where zeros denote
 corners in one tetrad and 1 denotes corners in the other; the last bit
 can be ignored, as it can easily be deduced from the other 7). This is
 slightly less space-efficient than computing the actual corner-separation
 coordinate as a value between `0` and `C(8,4)`, but it is much faster.
 
 We can now access the value in the `cocsepdata` table corresponding to
 our cube position. This table contains the value of the symmetry-reduced
 coordinate and the so-called *ttrep* (transformation to representative)
 necessary to compute the full **h0** coordinate. For convenience,
 this table also stores a preliminary pruning value based on the corner
 coordinate.  If this value is large enough, we may skip the computation of
 the **h0** coordinate and any further estimate. These three values can be
 stored in a single 32 bit integer, so that the table uses less than 12MB.
 
 ### Getting the pruning value
 
 Once we have computed the full h0 coordinate, we can access
 the correct entry in the full pruning table. We chose to
 use 2 bits per entry, as described by Rokicki in the [nxopt
 document](https://github.com/rokicki/cube20src/blob/master/nxopt.md).
 This means that every pruning table needs to have an associated *base
 value*, that determines the offset to be added to each entry (each entry
 can only be 0, 1, 2 or 3).  If the base value is `b`, a pruning value
 of 1, 2 or 3 can be used directly as a lower bound of b+1, b+2 and b+3
 respectively. However, a value of 0 could mean that the actual lower
 bound is anything between 0 and b, so we cannot take b as a lower bound.
 
 #### Fallback tables
 
 To be able to still use some sort of pruning value even when we get a
 0 read, we use a **fallback table**. Inspired by nxopt, this table is
 interleaved with the main table for cache efficiency: every 254 entries
 (508 bits), we store in 4 bits the minimum of the pruning values of these
 254 entries as a number from 0 to 16, without subtracting the table's
 base. Since a cache line is 512 bits long, we never get a cache miss
 when looking up the minimum value in the fallback table after a 0 read.
 Smaller lines (of 256 or 128 bits) have been tried, but they do not
 give any significant improvement over 512 bit lines.
 
 This trick provides the gratest performance benefits if the main pruning
 table is properly aligned. Unfortunately, as a design choice, the
 solver is implemented here as a library that defer all memory allocation
 business to the implementor. We do make sure that the table is properly
 aligned in the programs provided in this repository (for example, the
 rudimentary shell and the tools), but we do not enforce this in the
 main library code.
 
 #### Additional (fast) lookups
 
 Moreover, as an additional heuristic, in case of a 0 read we also look
 up another pruning value in a table that takes into account only the
 position of the edges. This table is small (around 1MB), so repeated
 accesses to it are not too slow. In practice, this gives a small speed up
 of around 5% for random scrambles. However, for small solvers (low values
 of h) and positions where corners are close to solved, this fallback
 table gives dramatic improvements (up to a factor of 1000x in some
 manual tests).
 
 More tables could be used to refine the fallback estimate, but each
 additional table leads to longer lookup times, especially if it is too
 large to fit in cache.  Previous versions of this implementation (up
 to commit 6c42463, or to version 0.2) also included the possibility of
 a table with 4 bits per entry, at least for the `h0` case. Such tables
 did not require a fallback lookup, but due to their large size they were
 not less efficient. Therefore, they have been removed.
 
 ### Estimation refinements
 
 After computing the pruning value, there are a number of different tricks
 that can be used to improve the estimation.
 
 #### Inverse estimate
 
 A cube position and its inverse will, in general, give a different
 pruning value. If the normal estimate is not enough to prune the branch,
 the inverse of the position can be computed and a new estimate can be
 obtained from there.
 
 #### Reducing the branching factor - tight bounds
 
 *This trick and the next are and adaptation of a similar technique
 introduced by Rokicki in nxopt.*
 
 We can take advantage of the fact that the **h0** (and **h11**) coordinate
 is invariant under the subgroup `<U2, D2, R2, L2, F2, B2>`. Suppose we
 are looking for a solution of length `m`, we are at depth `d` and the
 pruning value `p` *for the inverse position* is a strict bound, that is
 `d+p = m`.  In this situation, from the inverse point of view the solution
 cannot end with a 180° move: if this was the case, since these moves are
 contained in the target group for **h0** (or **h11**), it must be that
 we were in the target subgroup as well *before the last move*, i.e. we
 have found a solution for the **h0** coordinate of length `p-1`. But this
 is in contradiction with the fact that the inverse pruning value is `p`.
 
 From all of this, we conclude that trying any 180° move as next move is
 useless (because it would be the last move from the inverse position).
 We can therefore reduce the branching factor significantly and only try
 quarter-turn moves.
 
 #### Reducing the branching factor - switching to inverse
 
 We can expand on the previous trick by using a technique similar to NISS.
 
 Suppose that we have a strict bound for the *normal position*. As above,
 we can deduce that the last move cannot be a 180° move, but this does not
 tell us anything about the possibilities for the next move. However, if
 we *switch to the inverse scramble* we can take advantage of this fact
 as described above. For doing this, we need to replace the cube with its
 inverse, and keep track of the moves done from now on so that we can
 invert them at the end to construct the final solution.
 
 When this technique is used, it is also possible to avoid some table
 lookups: if the last move applied to the cube is a 180° move *on the
 inverse position*, then the coordinate on the normal position has not
 changed. Thus if we keep track of the last computed pruning value,
 we can reuse it and avoid an expensive table lookup.
 
 ### Pruning pipeline and prefetching
 
 To improve the memory access pattern and exploit
 [prefetching](https://en.wikipedia.org/wiki/Cache_prefetching)
 opportunities, we don't expand each neighbor one by one in the
 pruning phase. Instead, we employ a *pruning pipeline*, where we
 first compute the cube and inverse position for each neighbor, and
 then we proceed through a series of stages, at each of which we
 compute some prune off some nodes and we pre-compute the data
 necessary for the next step.
 
 This pipeline-based strategy gives an speedup of around 25%.
 Adding manual prefetching, we get an additional 20% speedup.
 
 ### Other optimizations
 
 The H48 solver uses various other optimizations.
 
 #### Avoiding inverse cube computation
 
 Computing the inverse of a cube position is expensive. We avoid doing
 that (for the inverse pruning value estimate) if we bring along both
 the normal and the inverse cube during the search, and we use *premoves*
 to apply moves to the inverse scramble.
 
 #### Multi-threading
 
 The solution search is parallelized *per-scramble*. Although when solving
 multiple positions it would be more efficient to solve them in parallel
 by dedicating a single thread to each of them, the current H48 solver
 optimizes for solving a single cube at the time.
 
 To do this, we first compute all cube positions that are 4 moves away
 from the scramble. This way we prepare up to 43254 *tasks* (depending
 on the symmetries of the starting position, which we take advantage of)
 which are equally split between the available threads.  Any solution
 encountered in this step is of course added to the list of solutions.
 
 #### Filtering out symmetric cases
 
 If the cube position we wish to solve is, say, mirrored along the RL
 plane, it is redundant to attempt starting a solution with both R and L',
 as they would lead to mirrored versions of the same solution. Therefore,
 we filter out some tasks by looking for symmetries before each of the
 first 4 moves. Then, unless we are looking for only one solution, we
 have to make sure to report back all symmetric variations of a solution,
 so we need to keep track of which transformations were available at each
 of these 4 points.
 
 This can reduce the number of starting positions drastically: for
 example, for the superflip our method produces 1038 starting tasks -
 a 41.67x improvement! Of course this optimization has no effect on any
 position that isn't 3 moves or fewer away from a symmetric position -
 which is the vast majority of them.
 
 #### Heuristically sorting tasks
 
 The tasks described in the previous paragraph (multi-threading) are
 initially searched in an arbitrary order. However, after searching
 at a sufficient depth, we have gathered some data that allows us to
 make some heuristical improvements: the tasks that leads to visiting
 more positions (or in other words, where we go over the estimated lower
 bounds less often), are more likely to yield the optimal solution. Thus
 we sort the tasks based on this, adjusting by a small factor due to the
 fact that sequences ending in U, R or F moves have more continuations
 than those ending in D, L or B moves - as we don't allow, for example,
 both U D and D U, but only the former.
 
 Preliminary benchmarks show a performance improvement of around 40%
 when searching a single solution. When searching for multiple optimal
 solutions the effects of this optimization will be less pronounced, and
 they are obviously inexistent when looking for *all* optimal solutions.
 
 ## Pruning table computation
 
 We first describe how the pruning table would be computed if we stored
 each pruning value fully, using 4 bits per entry. This algorithm used
 to be implemented (up to commit 6c42463), but has since been removed.
 
 ### Full tables for h0 and h11 (not implemented)
 
 If we are allowed to store the full pruning value for each entry,
 computing the pruning tables is not that difficult.  First, we set the
 value of the solved cube's coordinate to 0 and every other entry to 15,
 the largest 4-bit integer. Then we iteratively scan through the table and
 compute the coordinates at depth n+1 from those at depth n as follows:
 for each coordinate at depth n, we compute a valid representative for it,
 we apply each of the possible 18 moves to it, and we set their value to
 n+1. Once the table is filled or we have reached depth 15, we stop.
 
 Unfortunately what I described above is an oversimplification: one
 also must take into account the case where the corner coordinate is
 self-symmetric, and handle it accordingly by making sure that every
 symmetric variation of the same coordinate has its value set at the
 right iteration.
 
 When this algorithm reaches the last few iterations, the coordinates at
 depth n are many more than those whose depth is still unknown. Therefore
 it is convenient to look at the unset coordinate and check if they have
 any neighbor whose depth is i; if so, we can set such a coordinate to i+1.
 Thanks to [torchlight](https://github.com/torchlight) for suggesting
 this optimization.
 
 Finally, this algorithm can be easily parallelized by dividing the set
 of coordinates into separate sections, but one must take into account
 that a coordinate and its neighbors are usually not in the same section.
 
 ### 2 bit tables
 
 Computing the pruning tables for intermediate coordinates is not as
 simple.  The reason is that these coordinates are not invariant under
 the 48 transformations of the cube, but we are treating them as such.
 Take for example the case of **h10**. Each **h10** coordinate is
 represented by either of two **h11** coordinates: the one where the 11th
 edge is correctly oriented and the one where this edge is flipped (of
 course, the orientation of the 12th edge can be deduced from the parity
 of the orientation of the other 11). We can take any cube whose **h11**
 coordinate is one of these two to represent our **h10** coordinate, but
 the set of neighbors will change depending on which representative we
 picked. In other words: if I take one of the two cubes and make a move,
 the position reached is not necessarily obtainable by applying a move
 to the other cube. This means that the same algorithm that we described
 for the "real coordinate" case cannot be applied here.
 
 It is still possible to do a brute-force depth-first seach to fill
 these tables.  To make it reasonably fast, we apply a couple of
 optimizations. First of all, we restrict ourselves to computing a 2 bits
 table, so we can stop the search early.
 
 The second optimization we employ consists in pre-computing all possible
 **h11** coordinates at a fixed depth, and storing their depth in a hash
 map indexed by their **h11** coordinate value. This way we do not have
 to brute-force our way through from depth 0, and it make this method
 considerably faster. Unfortunately, since the number of coordinates
 at a certain depth increases exponentially with the depth, we are for
 now limited at storing the coordinates at depth 8.
 
 This method works for the "real" coordinates **h0** and **h11** as well.
 Moreover, in this case one can optimize it further by avoiding to repeat
 the search from a coordinate that has already been visited. Further
 optimization are possible for **h0** and **h11**, and we may implement
 them in the future.