Solving Full N × N × N Rubik ’ s Supercube Using Genetic Algorithm

The article presents an algorithm that uses an evolutionary approach to the problem of solving the Full Rubik N × N × N Supercube, i.e., the orientation of all cubies, including the internal ones, not only according to face colors but to the same orientation in 3D space. The problem is formally de ﬁ ned by the matrix representation using a ﬃ ne cubies transforms. The Full Supercube ’ s solving strategy uses a series of genetic algorithms that try to ﬁ nd a better cube con ﬁ guration than the current one. Once found, movements are made to change the current con ﬁ guration. This strategy is repeated until the cube is solved. The genetic algorithm limits the movements to the current cluster by solving the cube in stages, outwards from the center of the cube. The movements that solve the clusters are saved as macros and used to train and speed up the algorithm. The purpose of the presented algorithm is to minimize the solution time and not necessarily the number of moves.


Problem
The problem under consideration concerns solving new type of the N-segment Rubik's N × N × N supercube, which includes orientation of all internal cubies, using the genetic algorithm. This problem is much more general than solving traditional cubes due not only the need to correct orientation of the off-border cubies (and not only to the appropriate face color) but also to orient the inner cubies. To our knowledge, such structures, though being natural generalizations of Rubik's cube, were not yet investigated.
The traditional cube has three layers, it does not contain inner cubies, and the orientation of the middle cubies is usually not important. The first solving algorithms were developed in the early 1980s, and the most popular is Thistlethwaite's [1], based on the division of the task into 4 stages, in which the number of possible moves is gradually reduced. These restrictions apply to the possible angles of rotation for certain segments to only 180°. In the last step, all segments are subject to this limitation. For cubes larger than traditional ones, it is more important to limit the move-ments of the increasing number of segments, rather than types of rotation axis or values of rotation angles.
The problem of arranging a traditional cube with the use of the shortest possible sequence of movements is often solved by algorithms for iterative A * searching of the state tree, down to a specific depth (e.g., IDA * ). A typical example is Kociemba's algorithm [2], which also breaks the problem into substages but uses state tree search and heuristic functions. The most common approach, however, is to use special macromovement sequences (macros) that change the orientation of a small number of cubies, smaller than that resulting from a single movement [3]. Recent research has allowed to determine the minimum number of moves necessary to arrange N × N × N standard cubes (without internal cubies) from any configuration to OðN 2 /log NÞ, although the optimal algorithm is not yet known [4,5]. Already in 1994, Herdy and Patone solved the cube using a genetic algorithm that used macromovements as genetic operators [6]. Later studies tried to develop strategies without using expert knowledge [7,8] or to use already known heuristic algorithms [9,10].
The aim of the presented work was to develop the most effective and stable algorithm to solve the presented problem in the shortest possible time (number of GA), regardless of the number of movements made. With the problem posed in this way, the presented solution is based on the execution of a sequence of genetic algorithms improving the configuration of the cube in accordance with the adopted evaluation function and defined genetic operators. The solution found by the genetic algorithm is used by executing movements written in the chromosome of the best individual only if it improves the cube configuration; otherwise, the attempt is repeated. The solving algorithm is divided into stages, and the sequences that complete the stages are memorized and used for teaching and speeding up of the algorithm.

Cubies Orientation
Rubik's N × N × N cube has been defined in a coordinate system whose center coincides with geometric center of the cube. Each cubie of a cube is related to its local coordinate system defined in relation to the parent system of the entire cube. Movements of cube slices, in such defined system, result in additional rotations R of the cubies around their selected axis. Each cubie from the slice changes its origin position from p to p ′ according to the linear formula: After transformation, cubies swap their position and change orientation; therefore, configuration of the cube can be represented as a 3D array of cubies composite transformations. All cubies that can occupy the same position as a result of rotation, belong to the same set, called a cluster.
Positions p of the cubies do not have to be stored in the matrix, because they are dependent on their indices v in the 3D matrix. Since the indices are always nonnegative in the range ½0 : N − 1, the cubie coordinates p can be defined, e.g., by shifting their indices v by the indices of the cube center C, so that a shifted center will move to the origin of the cube coordinate system: where Positions p would be then origins of the cubies in the range ½−ðN − 1Þ/2, ðN − 1Þ/2 in each dimension, with an offset step of 1 from one another. The matrix of their transformation is then a combination of the translation to the center of the cube's system, the orthogonal rotation matrix, and inverse of this translation: The orientation and position of the cubies are thus clearly defined by the hitherto position of the cubie p and the rotation matrix R, which can be expressed as the composition of basic transformation matrices of the cube slice rota-tion in the planes of its system, around the X, Y, or Z axis by 90°. In column notation, these matrices are defined by Performing successive rotations in the same axis results in the combination of the same transformations, i.e., the exponentiation of these matrices. Obviously, their fourth powers are identity matrices (full angle rotations), hence, Since four orientations with the face in front are possible for each face of the cube (four multiples of the angle of 90°a round the vector perpendicular to the face), the total number of possible cubie orientations (as well as the maximum number of different positions of the cubie p′ or the maximum number of cubies in a cluster) is 24. Therefore, these orientations can be written shortly as a 5-bit number, but the most convenient way is to define them by using the Euler angles, presenting rotation R as a combination of elementary rotations around the cube axis in the order X, Y, and Z, respectively, by the angles ½α ′ , β ′ , γ ′ = π/2½α, β, γ, where α , β, γ = 0 : 3.
The determination of orientation by means of the Euler angles is not unequivocal. Different combinations of angles may specify the same orientations. For any rotation Q, a conjugation, like changes only rotation axis from axis to axis ′ = Q −1 axis. A real eigenvector axis ′ corresponds to the eigenvalue λ = 1 2 International Journal of Computer Games Technology and does not change its direction during rotation Q −1 R axis Q, and hence, it is an axis of this rotation: Using as a Q, e.g., 90°rotation about the X axis, changes rotation axis from the Z axis to the Y axis, because Therefore, the following equations are true: For α, β, and γ taking values in the range [0: 3], all possible cubie orientations can be obtained with just two movements; i.e., one of the angles can be considered to be equal zero: Which orientation definition can be used may be determined, e.g., by the position of the zero element in the rotation matrix (in bold in equation (11)). Due to the adopted order of transformations, rotations in the form R Y are enough to determine any orientation. We chose the first option. Using equation (10), it can be proved that for every possible value of β, the matrix R can be represented as a combination of just such rotations: Since the angles are multiples of 90°, it is easy to determine them from the signs of the sine and cosine functions ( Table 1, with φ = α, β, γ), which are elements of the rotation matrix ( Figure 1).
εis some small nonnegative value for taking account of float-point errors (smaller than cubieScale from Figure 2).
From the point of view of the movements performed, more favorable orientations are those that are composed of fewer moves. Any orientation can be defined by two or less moves. The orientations obtained as a result of a single movement must have one positive value in the cubie rotation matrix (equal to one in the absence of cubie scaling) on its diagonal. This number can also be determined directly from the values of γ, β, and α. Accurate orientation information is also important as it will allow evaluation function to recognize and remember cube configurations. The cubie's orientation vector, briefly referred to as its state, was encoded on 8 bits. Each of the consecutive 2 bits represents the number of turns (see Table 2): Orientation R 2 z R 2 X can be more preferably achieved with one move. In this case, it is worth to change its definition to R 2 Y . Extraction of the Euler angles from cubie composite affine transform R and definition of its orientation state are presented also in a form of pseudocode ( Figure 3).   Figure 2: Arrangement of cubies in the cube coordinate system. Table 2 presents possible cubie orientations in the form of its rotation matrix. For any cube configuration, this matrix can be recreated on the basis of a 3 × 3 submatrix of cubie's affine transformations-without taking into account its position (translation). Because each cubie is also scaled by positive valued cubieScale (Figure 2), only meaningful signs of matrix elements are presented.

International Journal of Computer Games Technology
The cube configuration can be saved as the orientation table of all cubies. Due to the small number of orientations of a single cubie (i.e., 24), they can be encoded as alphanumeric characters and the cube configuration as a string. This makes it easy to compare different cube configurations.
For a cubie p = ½x, y, z, the remaining cubies in the cluster are p′ cubies with coordinates Rp given in Table 2. Solv-ing the cube corresponds to the situation when for all cubies their rotation matrices are identity matrices-if the cubie is correctly oriented, it must be in the correct position.
Since the orientation of the inner cubies is important and the arrangement of the cubies begins from its center, all cubies are separated from each other by a certain constant spacing and scrambled cubies are presented as partially transparent ( Figure 2).  Table 2: Possible orientations of the cubie in the form of Euler angles (α′, β′, and γ′), matrix R, coded state, and cubie positions p′ after transformation by rotation R. state 00000000 01010000 01100000 01110000 01000100 01001100

Solution
International Journal of Computer Games Technology move slice, axis, angle ð Þ = 9slice + 3axis + angle: The gene values can then also be easily decoded into these 3 parameters, determining the movement.
The cube configuration can be represented as a block vector p, elements of which are the cubies rotation (orientation) matrices. For a single movement, permutation of the cubie p in the block vector can then be described as where the matrix M is the transformation of the rotation by the angle around the axis for cubies in the slice, and the identity matrix for the others: Since all cubies in the slice segment are rotated (permutated) in the cube vector, the cubie's transformation can be determined based on the M transformation combined with the transformation stored in the cube vector at the position from inversed transformation performed on that position: By saving the cube configuration as a vector, it is also possible to define a transformation that performs any permutation of the elements of this vector (cubies in a cube), composite with the M transformation, and present the movement of the cubie slice as a matrix operation: Matrix ℝ has the size N 3 × N 3 , and its elements are the M matrices of the cubies. There is one nonzero M transformation at the position specified by undoing the M transfor-mation at the original cubie position and converting it to a linear index in each row and column. Knowing the current cubies orientations, it is possible to construct such a matrix ℝ, which is a composite of all previous cubie's transformations and represents cube configuration (examples for N = 3 are presented in Figure 4). For example, for a cube with size N = 2 and unfolding to the vector p = ½p XYZ , the single rotations of the slice = 0 around the Z and Y axes are represented by matrices (dots represent zero matrices): , and commutator It follows that only 4 cubes changed their transformations, but they swapped their places p 000 − p 010 and p 100 − p 101 and are transformed by a single move. Executing a

5
International Journal of Computer Games Technology double commutator will result in the absence of a permutation of the positions of these cubies with composite transformation from 2 movements: Thanks to this, the movements of commutators can be used to position the corner cubies of the cube (see also Figure 4).
The cube configuration changes through consecutive movements: Finding the movements that reverse these transformations, that is, solving the cube, is related to the decomposition of the matrix ℝ T into the sequence of the cube movements: Each of these transformations is orthogonal to a given segment. As in the one-side Jacobi transform used to derive the SVD of a matrix by means of a Q T RQ transformation to a similar matrix (maintaining the same eigenvalues), such a sequence in the case of a cube transformation, called rotation conjugate of R, also changes only a limited number of elements. Presented algorithm will also use a strategy similar to Gaussian elimination or "chasing zero" in QR iterations, by solving the cube, cubie by cubie in a specific order, so as not to cause deterioration of the existing configuration.
The number of possible configurations of a cube with N segments can be estimated, by assuming that each cubie can be in any orientation, on This number therefore increases exponentially very quickly, and the chromosome evaluation function must take into account not only the orientation of the N 3 cubies themselves but also their position, changed as a result of rotation.
Due to the possibility of repeating the same cubie position coordinates during rotation, there are 6 types of 6 International Journal of Computer Games Technology clusters. They differ in the number of segments that can affect the position of cubies of a given cluster type. For each axis and selected angle, the index of the segment influencing the cubie in the cluster can take the value of the cubie coordinate p ' after transformation, i.e., the cubie coordinate x, y , z of p or the opposite number (Table 2), which gives a maximum of 6 options: Apart from the possibility of changing the axis and the angle of rotation for each segment, the type of the C k cluster can be determined precisely by the number of k = 1::6 segments affecting it. Figure 5 shows all possible cubies positions from given C k cluster as red squares in a transparent cube. Vertical lines mark segments for chosen axis. If we try to arrange a cubie or cubies only from a given cluster, then the values of genes in the chromosome can be restricted to movements that affect these cubies. The number of such movements for any cube size, taking into account the possibility of choosing one of the three axes and one of the three angles, is 9k (max. 54).
Clusters C 2 consist of 8 corner cubies for cube or for inner cubes. For odd-sized cubes, the C 1 cluster constitutes one, middle cube, and the C 3 clusters consist of 6 cubies located in the middle of the cube segments (except for the middle segments) or 12 corner cubies for the middle segments of the cube or inner cubes. The remaining types of clusters contain maximal number of cubies, i.e., 24 (presented in Figure 6). Even N size cubes do not contain odd cluster types.
A set of cubies in a cluster can be created on the basis of a selected cubie p by applying to it any possible rotation R and obtaining positions (or orientations) of cubies p ′ , which belong to the same cluster (Table 2 and Figure 7).
Based on restricted slice indexes and any possible value of axis and angle, list of restricted moves is created. It serves as a gene pool for creating initial population of chromosomes in genetic algorithm. The chosen genetic operators will not violate this restrictions, and whole algorithm will run on restricted list of acceptable moves. Implementation example of function retrieving restricted list of moves is presented in Figure 8.
Because in the definition of movement, as the size of a cube N increases, only the number of segments increases; the number of possible moves for a cube should be limited by blocking the movements of specific slices. Additionally, when the slice is rotated, cubies do not change their distance to the cube center. Movements of the slices in planes that are more distant from the center than r do not affect the position of cubies in planes less distant. Thus, these cubies can be arranged only by means of slice movements that are distant from the center by less than r. The adopted strategy solves cube by orienting cubies one by one, preserving orientation of the cubies previously arranged, whereby selected (active) unoriented cubie is always that, which is closest to the middle of the inner cube of the smallest possible size. This way, the algorithm will try to arrange the inner cubies   International Journal of Computer Games Technology from inside to the outside of the cube. The algorithm tries to correctly orientate active cubie by limiting the movements to the movements influencing its cluster. By arranging the cube from inside, one can gradually limit the movements of the interior slices. The last clusters of the cube (and also inner cubes) are always the clusters containing the corner cubies-most distant from the center. The movements of these clusters are restricted to only movements of the 6 outermost walls of a given cube. The consecutive stages of solving 5 × 5 × 5 supercube are presented in Figure 9.
The whole process of solving the 3 × 3 × 3 cube, in form of ever changing configuration matrix, is presented in Figure 10. White squares represent zero matrices, and black ones represent identity matrices. Identity block matrix is a configuration of a solved cube.

Fitness Function.
A chromosome containing a sequence of cube movements may contain movements that can be combined into a single movement or that are completely self-canceling. Prior to the evaluation of the chromosome, a correction procedure can be performed that combines the movements of the same slice of the cube into one movement if they are not separated by movements in different planes, or to remove both genes when movements are cancelling each other.
The fitness function evaluates the cube's alignment status for each subsequence of movements in the chromosome starting with the first move. The first movement on a chromosome should change the orientation of the active cubie, changing the number of the movement's slice to the segment containing the active cubie, taking into account the axis of rotation specified in the gene. This will increase the probability of finding the active cubie orienting sequence.
First, a copy of the cube is created, which will be subject to movements written in the chromosome. The chromosome genes are processed sequentially. The gene is decoded to the movement and performed on a copy of the cube, which is then evaluated. The sequence of movements is successively supplemented with successive movements defined in the genes of the chromosome, creating an ever longer substring all the way to the entire chromosome. The length of the substring with the best score is remembered, and its score becomes the overall chromosome score. In this way, the possibility of finding solutions with a different number of movements using chromosomes of a fixed length has been effectively implemented.
The cube is assessed on the basis of the sum of the assessments of individual cubies belonging to the selected set of cubies. This collection consists of two groups. Group A consists of all cubies belonging to solved clusters. Position of cubies belonging to this group cannot be changed by the algorithm; therefore, they gain a weight greater than the sum of the assessment of all cubies from group B, which are cubies belonging to the current active cubie cluster. The value of 0 is the best fitness of the cube arrangement and individual cubies. The fitness of the X = A or B group is the sum of the scores of nonoriented cubies, each in the range [0, 2) and defined as where X i denotes orientation of the i-th cubie, encoded as its state (see equation (13) and Table 2) in A or B group. Cube fitness error is a combination of both group fitnesses: where the maximum value of group B's fitness is the double number of cubies in this group (according to equation (25)).
By taking into account the exact orientations of the cubies, the fitness function causes the genetic algorithm to prefer placing cubies in a certain order, always bringing the After arranging cubies from clusters in groups A and B, the value of the cube fitness function drops to zero. A new active cubie is designated. If the cube is not yet solved, then such a cubie exists and it is located at a distance from the center of the cube greater than the cubies in groups A and B. The group B is then joined with A, and the cubies from the active cubie cluster form a new group B. This ends the cube building stage, and the solution that led to it, if not a sequence derived from a saved solution, is stored in the solution file and can be used to speed up the algorithm. The cube's fitness is set to its maximum value, and another genetic algorithm is executed.

Selection, Crossover, and
Mutation. An important step in the genetic algorithm is selection. In the case of a Rubik's cube, the chromosomes will often differ in genes at distant positions in the chromosome, which will not affect the value of its fitness function. When selecting winners, a kind of ranking selection can be used, which omits individuals with  International Journal of Computer Games Technology the value of the fitness function the same as the already selected winner. The advantage of this approach is also low computational complexity.
The most important thing in the algorithm, however, is the use of efficient genetic operators. While the crossover may be a simple one-point recombination of chromosomes, the mutation should allow finding a solution for a cube in a configuration of only a few cubies with wrong orientation, i.e., macros. Such R′ sequences are in the form of a commutator or a conjugation in which the matrix Q can be composed of several movements around different axes: By limiting the movements to the movements that change the active cubie's cluster, it is much easier to find the sequence improving the cube's configuration in accordance with the fitness function. If the R axis performs rotation on a slice more distant from the center of the cube than the inner slices of the already oriented cube, then it is an identity transformation for the cubies of the inner cube. This means   2  31  99  130  364  GA2  2  21  88  190  459  GA3  4  17  85  175  405  GA4  3  17  44  280  727  GA5  4  16  46  144  493  GA6  3  18  95  134  545  GA7  2  16  84  256  836  GA8  3  13  69  225  776  GA9  2  14  66  184  507  GA10  3  14  112  166  493   10 International Journal of Computer Games Technology that also the entire R ′ transformation is an identity transformation for these cubies and does not change their arrangement.
The proposed mutation procedure performs the following actions: (1) Random selection of an idx less than half-length of the chromosome (2) Genes are processed from the start of the chromosome to idx − 1, and inverse transformations to them are coded in genes from index 2•idx down to idx + 1 Due to the high efficiency of finding macro sequences, the mutation rate in the algorithm can be extremely high, within 100% of the population size.

3.4.
Learning. The sequences found by the genetic algorithm that end with the arrangement of all cubies from the assessed set cause the change of position of only a few, k ≪ N 3 cubies, not changing orientations of the rest. These are the key sequences that can be used not only in the last stage of solving the cube. The algorithm can remember such sequences, implementing the ability to learn and arrange the cube faster and faster.
The memorized sequences are used each time a better cube configuration is found to further improve it. The movements of these sequences are restricted to a certain cluster (maximum 6 segments) and can be mapped to active cubie cluster movements by replacing segments with active cluster segments. This enables, among others, using the saved movements from smaller cubes, not only to orient the inner cubes but for all clusters of the same type in the cube.
Movements that, after composing, can be represented in the form of an orthogonal matrix R can be additionally complemented by a conjugation to the form Q −1 RQ, in which case the matrix Q is a single rotation, selected optimally among all possible movements of the active cluster. This will enable a much wider use of memorized movements, taking into account similar patterns resulting from the cube symmetry.

Results
The way the algorithm works can be summarized in four main points:  as an active cubie that, which is closest to its middle and contained in the inner cube of the smallest possible size (2) Gene values in chromosomes (movements) are limited to segments that influence the current cubie cluster, i.e., the current cluster (3) Only cubies belonging to the current cluster are assessed. Cubies from oriented clusters cannot change orientation (4) The algorithm has the ability to learn macros, which are the sequence of movements that complete orientation of the current cluster. Macros are mapped to the movements of the current cluster by changing segment indexes to segment indexes of the current cluster The results of the algorithm's operation on full cubes with sizes from 2 to 6 are presented. Each cube was scrambled with a 100 random moves before the algorithm was started.
The parameters of the algorithm assumed the length of the chromosome equal to 27, and the number of selection winners was equal to 10% of the population, the population itself at the level of 100 individuals, and the number of generations at the level of 50 iterations, due to the attempt to find only the cube configuration better than the one found so far. For the duration of the experiment, the ability of the algorithm to learn was suspended. The algorithm has had already known and used 100 or 200 macro sequences. Tables 3 and 4 ( Figure 11) and Tables 5 and 6 ( Figure 12) show the number of genetic algorithms performed and the time it took to complete the cube in 10 consecutive experiments.

12
International Journal of Computer Games Technology Tables 7 and 8 and Figure 13 show the number of moves, which led to a solution. Some of them could stand from learned macros.

Conclusions
The structure of the Full Rubik's Supercube differs from the classic supercube [11]-and the existing algorithms (reduction and cage) and macromovements, which are used to solve it, may cause changes in the orientation of the inner cubies. This makes the puzzle much more difficult, and comparing the times or the number of moves needed to solve a full and classic supercube is not directly possible for N > 3. For N = 3, a full supercube differs from a supercube only in the need to orient one inner cubie. Comparison of the results, for example with the algorithms presented by El-Sourani et al. [9], shows a greater number of moves necessary to solve a cube in our algorithm, but much shorter times to find a solution (much more important when solving larger cubes). To our knowledge, the Full Rubik's Supercube, from the point of view of its solving, has not yet been analyzed. However, it is certainly possible to determine the strategy and sequences of macromoves for full supercube as well. This task is well fulfilled by the methodology presented in the article, which uses a sequence of genetic algorithms in combination with a simple technique of learning detected macros.
Presented method deals with a generic problem of solving Full Rubik's Supercube of any size, taking into account orientation of the off-border and internal cubies. The cube is represented as a model of the composite transformations of the block rotation matrices around X, Y, and Z axes. The resulting cube configurations are therefore always correct, and it is not necessary to color map the faces of cubies and select correct permutations to determine their orientation. The article presents also a way of formally writing the cube configuration as a combination of moves (transformations) in the form of a block matrix. This type of notation can be helpful especially for the analysis of generalizations of the Rubik's cube in multidimensional spaces. In n -dimensional space, the orientation of the cubies is given by n 2 ! Euler angles. As the number of dimensions increases, the number of possible movements therefore also increases, due to the increasing number of possible planes of rotation. Future research will focus on the greater use of cube symmetry and the problem of solving a multidimensional cube. As we mentioned, the problem of solving the cube can be formulated as a decomposition of configuration block rotation matrix into orthogonal move block matrices (see also equation (23)). Such analysis can be treated as a generalization of the SVD decomposition for higher dimensions [12], and similarities of the presented algorithm to this algorithms are therefore highlighted. Solving strategy, which arranges cubies in the current cluster, one by one, allows for strong restriction of the possible moves (gene values in GA), which depends on the type of cluster, but is indifferent to the size of the cube (number of segments). The method also does not use prior expert knowledge of the existing macromoves until it finds them on its own.
The results of all experiments were given for algorithms with a small number of learned macromoves (100 and 200) and the learning module turned off, in order to demonstrate the ability of the system to learn quickly and to properly compare the times of solving cubes of different sizes.
The project available in the repository was written in C# language with the use of OpenGL functions and without the use of other programming libraries. It has the characteristics of a game, because it is very satisfying to watch how the algorithm learns to solve a full cube faster and faster. It can be easily supplemented with a learning module for manual supercube solving, based on commutator sequences learned by the system. Cubies with the wrong orientations can be rendered semitransparent, revealing internal structure of a Full Supercube.
Rubik's cube can fascinate, but it is not only a toy; it is also used in data encryption [13][14][15][16], search algorithms in 3D structures, or solving mechanical problems [17][18][19][20][21][22]. It is an inspiration to build scientific analogies in the theory of groups, permutations, cycles [23], but also in 3D graphics. Like fractals, multidimensional regular polyhedra, or block tensors, it is a gateway to new dimensions and discovering rules and beauty in mathematics.

Disclosure
A preprint has previously been published [24].