Lispで調べるルービックキューブ

この記事は [Lisp Advent Calendar](https://qiita.com/advent-calendar/2019/lisp) 2019の18日目の記事です。


---

# TL;DR
* 前半は、ルービックキューブ(ただし、$ 2 \times 2 \times 2 $)を解く(6面が揃える)仕組みを説明します。
* 後半で、解くために必要な手順を調べるcommon lispのプログラムを説明します。
* なお、この記事は[第9回関西Lispユーザ会](https://kansai-lisp-useres.connpass.com/event/113242/)の発表内容に加筆修正したものになります。



# ルービックキューブ
ご存知の方も多いと思いますが、ルービックキューブは立体パズルの一つです。立体パズルの中でも老舗に近いものではないかと思います。通常ルービックキューブと言えば、$ 3 \times 3 \times 3 $ ですが、ここでは簡便のため$ 2 \times 2 \times 2 $ を扱います。以下、ルービックキューブや数学的記述などは[群論の味わい](https://www.kyoritsu-pub.co.jp/bookdetail/9784320019416)を元ネタに使っています。

## 基本的な考え
達人は10秒もあれば、解く(6面を揃える)ことができますが、その様子を見ていても、何をしているのか判らない人が多いでしょう。そのため、ルービックキューブを解くことは大変難しいことのように思えますが、単に解く(時間をかけてもよい)のであれば、それほど難しいことではありません。 

$ 2 \times 2 \times 2 $のルービックキューブは三色の小方体八個から構成されています。八個の小方体が位置と向きがどう変わっていくかを考えます。この際、ルービックキューブの向きを固定して考えます。次の図では、オレンジの面を前、白の面を上になるようにルービックキューブの向きにしたものです。


![cube-move-photo.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/323744/d829f68f-e7a4-5539-0123-c8c88e53b442.png)

この状態で、前面、後面、右面、左面、上面、下面を時計回りに90度回す手順をそれぞれF、B、R、L、U、Dとします。90度回す手順の組み合せを表わすには$\circ$を使います。例えば、右面を90度回した後に上面を90度回す手順は$R \circ U$となります(結果は次の図)。


![RU.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/323744/f37e41ca-3ede-0b68-ff29-2e353c6d04d6.png)

ここでいくつかの決め事をしておきます。

* 全く何も操作しないことを$E$とします(たし算での0、かけ算での1の相当)。
* 180度回すことは、90度回すことを2回実施することなので、$R \circ R$ (右面を180度の場合)ですが、これを$R^2$と表わします。
* 反時計回りに90度回すことは270度回すことになりますが、$R^{-1}$と表わします[^1]。

これらより、$$R \circ R^{-1} = E$$ や $$R^{4} = E$$となります。

----

さて、唐突ですが、6面が揃った状態で
$$R \circ U \circ R^{-1} \circ U^{-1} \circ F^{-1} \circ U^{-1} \circ F$$
という手順を実施すると、次の図になります。

![swap.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/323744/5f9161d0-662c-31c9-964a-ade4bd2a09b0.png)

6面全てが見えているのでないので、判りにくいとは思いますが、下の4つの小方体は最初の状態と同じです。上の4つ小方体は、右側の二つは入れ替わっています(手前のは向きは変わっています)。左側の手前は変化せず、奥の小方体は向きは変わっていますが、位置は変わっていません(この図だけでは判りにくいとは思いますが)。このように、二つの小方体の位置を交換する手順が存在します。また、全て小方体の位置を動かさず向きのみを変える手順も存在します。

二つの小方体の位置を交換する手順を繰返していくことで、全ての小方体を元の位置に戻すことができます。その後、向きを変える手順を繰返していくと最終的にルービップキューブを解くことができまます。すなわち、ルービックキュープを解くためには、

* 位置を変える手順
* 向きを変える手順

の二つだけを知っていればよいのです。たった二つです。簡単でしょう。

このような手順を求めるのが今回の目標です。

## 位置と向き
位置の変化を示すために、各小方体に番号を与えます。向きの変化を示すために、上面の下面に印を付けておきます。この印を基準参照印と呼びます。次の図では白(上面)と黄(下面)にシールを付けています。


![position.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/323744/cc3d0a10-93f3-2d9e-9074-4fb08fe09667.png)

位置については、番号がどう移動したかにより判断します。上図から、$R$を実行すると、


* $ 1 \to 3$ 
* $ 2 \to 1$
* $ 3 \to 4$
* $ 4 \to 2$

と移動します。これは数学では置換と呼び、次のように表現します。

```math
    \left(
    \begin{array}{cccccccc}
    1 & 2 & 3 & 4 & 5 & 6 & 7 & 8 \\\\
    3 & 1 & 4 & 2 & 5 & 6 & 7 & 8 
    \end{array}
    \right)
```

8個の要素の置換全体を「8次の置換群」と呼び、$S_8$と記します。また、左側の小方体($5,6,7,8$)は位置が変わっていないことも判ります(当然ですが)。



向きについては、移動先において基準参照印から、(時計回りで)120度の何倍分ひねられているかを考えます[^2]。

* $1$は元の$3$の基準参照印から240度ひねられているので2
* $2$は元の$1$の基準参照印から120度ひねられているので1
* 左側の位置が変わっていない小方体は、向きも変化していないので0

となります。8個の小方体のひねりをまとめて次のように表わします。

$$(2,1,1,2,0,0,0,0)$$
これの全体は3を法とする加法群の8次のベクトルとなり、$C_3^8$と記します。

# 手順の演算
以下の説明は非常に天下り的です(証明はしないということ、証明については[群論の味わい](https://www.kyoritsu-pub.co.jp/bookdetail/9784320019416) 補題 11.1.3参照)。

$\mathit{g},\mathit{h}$をルービックキューブの手順とする。手順に対応する置換を$\rho$ 、向きを$\mathit{v}$とする。手順$\mathit{g}$を実行後に手順$\mathit{h}$を実行する手順は、$g \circ h$で表わされ、

* 小方体の並び換えは $\rho(g) \circ \rho(h) $ (置換の積)
* 小方体の向きは $\vec{v}(g) + \rho(g)^{-1}(\vec{v}(h))$

となります。


$R \circ F$を例に具体的に計算します。

![RF.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/323744/40a0154b-c701-90e3-ca5f-cb0f2b89db5d.png)

位置の計算
![F.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/323744/650fe0cc-f6a7-1091-e697-1d6ba34549d3.png)

```math
\rho(R) = \left(
    \begin{array}{cccccccc}
    1 & 2 & 3 & 4 & 5 & 6 & 7 & 8 \\
    3 & 1 & 4 & 2 & 5 & 6 & 7 & 8 
    \end{array}
    \right)
```

```math
\rho(F) = \left(
    \begin{array}{cccccccc}
    1 & 2 & 3 & 4 & 5 & 6 & 7 & 8 \\
    8 & 2 & 1 & 4 & 3 & 6 & 7 & 5 
    \end{array}
    \right)
```

```math
\rho(R \circ F) = \rho(R) \circ \rho(F) = \left(
    \begin{array}{cccccccc}
    1 & 2 & 3 & 4 & 5 & 6 & 7 & 8 \\
    1 & 8 & 4 & 2 & 3 & 6 & 7 & 5 
    \end{array}
    \right)
```


向きの計算
![R-inverse.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/323744/23bdcfe4-ace4-574b-b075-666741f15d0a.png)


```math
\rho(R)^{-1} = \rho(R^{-1}) = \left(
    \begin{array}{cccccccc}
    1 & 2 & 3 & 4 & 5 & 6 & 7 & 8 \\
    2 & 4 & 1 & 3 & 5 & 6 & 7 & 8 
    \end{array}
    \right)
```

```math
    \vec{v}(F) = (1,0,2,0,1,0,0,2)
```

```math
\begin{eqnarray}    
\rho(R)^{-1}(\vec{v}(F)) = \left(
    \begin{array}{cccccccc}
    1 & 2 & 3 & 4 & 5 & 6 & 7 & 8 \\
    2 & 4 & 1 & 3 & 5 & 6 & 7 & 8 
    \end{array}
    \right)
    ((1,0,2,0,1,0,0,2)) \\
= (2,1,0,0,1,0,0,2)
\end{eqnarray}    
```

```math
\vec{v}(R) = (2,1,1,2,0,0,0,0)
``` 

```math
\begin{eqnarray} 
    \vec{v}(R) + \rho(R)^{-1}(\vec{v}(F)) = (2,1,1,2,0,0,0,0) + (2,1,0,0,1,0,0,\
2) \\ 
    = (4,2,1,2,1,0,0,2) \\ 
    = (1,2,1,2,1,0,0,2)
\end{eqnarray}
``` 

# Lispによる実装
ここからLisp Advent Calendarらしくなります。

ある手順を実行すると、位置($S_8$)と向き($C_3^8$)が決まるので、ルービックキューブの手順を表現するために、この二つを要素と持てばよいことが判ります。CLOSを使って、クラス`rubik2`とします。


```lisp
(defclass rubik2 ()
  ((vertex
    :initarg :vertex
    :documentation "頂点の位置(置換)")
   (v-twist
    :initarg :v-twist
    :documentation "頂点の向き(ベクトル)"))
  (:documentation "ルービップキューブの手順(2*2*2)"))  

(defmethod print-object ((obj rubik2) stream)
  (format stream "#vertex:~a,twist :~a"
          (slot-value obj 'vertex)
          (slot-value obj 'v-twist)))
```
`vertex`は置換です。置換もクラスにします。

```lisp
(defclass permutation () 
  ((row
    :initarg :row))
  (:documentation "置換の上段が1 2 ...と固定した場合の下段とする"))

(defmethod print-object ((obj permutation) stream)
  (format stream "~a" (slot-value obj 'row)))
```

置換や手順の演算($\circ$)をgenericなメソッドで実装します。名前は`*_`とします。$\rho(g)^{-1}(\vec{v}(h))$を$\rho(g)^{-1} \circ \vec{v}(h)$と考えて、`*_`とします。また、演算の結果は新しいインスタンスを作成することにします。

```lisp
(defmethod *_ ((x permutation) (y permutation))
  (let ((p1 (slot-value x 'row))
        (p2 (slot-value y 'row)))
    (make-instance 'permutation
                   :row (loop for l in p1 collect (nth (1- l) p2)))))

(defmethod *_ ((x permutation) lst)
  "リスト lstに 置換 xを適用する"
  (with-slots (row) x
    (loop repeat (length row)
        for i from 1
          collect (nth (position i row) lst))))

(defmethod *_ ((x rubik2) (y rubik2))
  (let ((move1 (slot-value x 'vertex))
        (move2 (slot-value y 'vertex))
        (twist1 (slot-value x 'v-twist))
        (twist2 (slot-value y 'v-twist)))
    (make-instance 'rubik2
                   :vertex (*_ move1 move2)
                   :v-twist (mod3 twist1
                                  (*_ (inverse move1) twist2)))))
```

$R^{-1}$などの逆演算も`inverse`という名前のgenericなメソッドで実装します。

```lisp
(defmethod inverse ((x permutation))
  (with-slots (row) x
    (make-instance 'permutation
                   :row (loop repeat (length row)                               
                              for i from 1
                              collect (1+ (position i row))))))

(defmethod inverse ((obj rubik2))
  (with-slots (vertex v-twist) obj
    (make-instance 'rubik2
                   :vertex (inverse vertex)
                   :v-twist v-twist)))
```

手順の逆演算では、向きは変わりません。$U$と$D$について考えると、90度単位で、どれだけ(max3ですが)回しても基準参照印の付け方から、向きが変わることはありません。それ以外の面についても、180度回しても向きが変わることはありません。$R^{-1}$は$R$を180度回したものになりますから、手順の逆演算では向きは変わらないことが判ります。

$E,R,U,F$を事前に定義します[^3]。

```lisp
(defvar E (make-instance 'rubik2
                         :vertex (make-instance 'permutation :row '(1 2 3 4 5 6 7 8))                                 
                         :v-twist '(0 0 0 0 0 0 0 0)))

(defvar R (make-instance 'rubik2
                         :vertex (make-instance 'permutation :row '(3 1 4 2 5 6 7 8))
                         :v-twist '(2 1 1 2 0 0 0 0)))

(defvar U (make-instance 'rubik2
                         :vertex (make-instance 'permutation :row '(1 2 5 3 6 4 7 8))
                         :v-twist '(0 0 0 0 0 0 0 0)))

(defvar F (make-instance 'rubik2
                         :vertex (make-instance 'permutation :row '(8 2 1 4 3 6 7 5))
                         :v-twist '(1 0 2 0 1 0 0 2)))
```
$R \circ U \circ R^{-1} \circ U^{-1} \circ F^{-1} \circ U^{-1} \circ F$を計算してみます。

```lisp
CL-USER> (reduce #'*_ (list R U (inverse R) (inverse U) (inverse F)
                           (inverse U) F))
#vertex:(1 2 4 3 5 6 7 8),twist :(0 0 0 1 0 2 0 0)
```
3と4が入れ替わり、3は(4の位置で)120度、6は(6の位置で)240度、向きが変わっています。合っているようです。

# 手順を調べる
手数を決め、その総手順について調べるという単純な方法を用います[^4]。

## 数え上げる
`(R 1)`は$R$、 `(R 2)`は$R^2$、`(R 3)`は$R^{-1}$とする。これらを1手順とみなして、数え上げる関数`enumerate`を考える。 `(R 1) (U 1) (F 1)`の2手順の結果は次のようになる。 

```lisp
CL-USER> (loop for k being the hash-key in (enumerate '((R 1) (U 1) (F 1)) 2)
               do
                  (format t "~a~%" k))
((R 2))
((R 1) (U 1))
((R 1) (F 1))
((U 1) (R 1))
((U 2))
((U 1) (F 1))
((F 1) (R 1))
((F 1) (U 1))
((F 2))
NIL
```

計算途中で同一手順のダブリをなくすために、手順をハッシュ表のキーとして扱っている[^5]。このため`enumerate`の結果もハッシュ表。

```lisp
(defun enumerate (items num)
  "itemsの要素同士のnum個の組合せをhash表で返す"
  (flet ((init-hash ()
           (loop with dp-hash = (make-hash-table :test #'equal);
                 for k in items
                 do (setf (gethash (list k) dp-hash) nil)
                 finally (return dp-hash))))
    (if items
        (loop for i from 1 to num
              for result = (init-hash) then (direct-product-hash result items)
              finally (return result)))))

(defun direct-product-hash (prev items)
  "prev(hash表)のキーにitemsをつなぐことで直積を計算する。結果はhash表"
  (if (null items)
      prev
      (loop for x being the hash-key in prev
            with result = (make-hash-table :test #'equal)
            do
               (loop for y in items
                     for head = (butlast x)
                     for tail = (shrink (car (last x)) y)
                     if (or head tail)
                     do (setf (gethash `(,@head ,@tail) result) nil))
            finally (return result))))               

(defun shrink (x y)
  "xとyが同じ面であれば、回転を足す。"
  (if (eq (side x) (side y))
      (let ((rot (mod4 (rot x) (rot y))))
        (unless (zerop rot)
          (list (list (side x) rot))))
      (list x y)))

(defun side (x)
  "ex. xが(R 3)ならRを返す"
  (first x))

(defun rot (x)
  "ex. xが(R 3)なら3を返す"
  (second x))
```

`mod4`の定義は後で示します。

## インスタンスに変換
数え上げの結果として、`((R 3) (U 1) (R 2) (F 1))`のようなリストですが、
手順を調べるにはインスタンスに変換させることが必要。

```lisp
(defun repeat (symbol n)
  (make-list n :initial-element symbol))

(defun pairs->move-list (pair-list)
  (loop for pair in pair-list
        append (repeat (first pair) (second pair))))
```

実行例

```lisp
CL-USER> (pairs->move-list '((R 3) (U 1) (R 2) (F 1)))
(R R R U R R F)
```

数え上げて、インスタンスに変換するためのマクロ。


```lisp
(defmacro moves->obj ((&rest move-list) n)
  (let ((seed (loop for i from 1 to 3
                    append (mapcar #'(lambda (x) (list x i))
                                   move-list))))
    `(flet ((sym->obj (sym)
              (cond ,@(loop for m in move-list
                          collect `((eq sym ',m) ,m)))))
       (loop for pairs being the hash-key in (enumerate ',seed ,n)
             collect (list pairs 
                           (loop for m in (pairs->move-list pairs)
                                 collect (sym->obj m)))))))
```

マクロ展開結果。

```lisp
CL-USER> (macroexpand-1 '(moves->obj (R U F) 7))
(FLET ((SYM->OBJ (SYM)
         (COND ((EQ SYM 'R) R) ((EQ SYM 'U) U) ((EQ SYM 'F) F))))
  (LOOP FOR PAIRS BEING THE HASH-KEY IN (ENUMERATE
                                         '((R 1) (U 1) (F 1) (R 2) (U 2)
                                           (F 2) (R 3) (U 3) (F 3))
                                         7)
        COLLECT (LIST PAIRS
                      (LOOP FOR M IN (PAIRS->MOVE-LIST PAIRS)
                            COLLECT (SYM->OBJ M)))))
T
```

## 7手順を調べる
先に実行例を示します。

```lisp
CL-USER> (time
          (length
           (setq *move7* (moves->obj (R U F) 7))))
Evaluation took:
  323.743 seconds of real time
  323.473099 seconds of total run time (322.747987 user, 0.725112 system)
  [ Run times consist of 3.051 seconds GC time, and 320.423 seconds non-GC time. ]
  99.92% CPU
  466,191,771,282 processor cycles
  585,900,976 bytes consed
  
503883
CL-USER> (time (length (setq *fix7*
                             (loop for x in *move7*
                                   for move = (reduce #'*_ (second x))
                                   if (and (fix-p move '(1 2 7 8))
                                           (not (=_ move E))
                                           (not (=_ move U))
                                           (not (=_ move (*_ U U)))
                                           (not (=_ move (inverse U))))
                                   collect x))))
Evaluation took:
  58.583 seconds of real time
  58.584552 seconds of total run time (58.412399 user, 0.172153 system)
  [ Run times consist of 1.153 seconds GC time, and 57.432 seconds non-GC time. ]
  100.00% CPU
  24 lambdas converted
  84,359,374,794 processor cycles
  4,470,724,848 bytes consed
  
220
CL-USER> (loop for (pairs o-list) in *fix7*
               for obj = (reduce #'*_ o-list)
               if (=_ obj '(1 2 4 3 5 6 7 8))
               do
                  (format t "~a~% ~a~%" (pairs->formula pairs) obj))
R*U*R^(-1)*U^(-1)*F^(-1)*U^(-1)*F
 #vertex:(1 2 4 3 5 6 7 8),twist :(0 0 0 1 0 2 0 0)
F*R^(-1)*F^(-1)*U^(-1)*R^(-1)*U*R
 #vertex:(1 2 4 3 5 6 7 8),twist :(0 0 0 1 2 0 0 0)
R^(-1)*U^(-1)*R*U*F*R*F^(-1)
 #vertex:(1 2 4 3 5 6 7 8),twist :(0 0 2 0 1 0 0 0)
F^(-1)*U*F*U*R*U^(-1)*R^(-1)
 #vertex:(1 2 4 3 5 6 7 8),twist :(0 0 2 0 0 1 0 0)
NIL
CL-USER> (loop for (pair o-list) in *fix7*
               for obj = (reduce #'*_ o-list)
               if (not (find 2 (slot-value obj 'v-twist)))
               do
                  (format t "~a~% ~a~%"
                          (pairs->formula pair) obj))
R*U*R^2*U^(-1)*R^2*U*R
 #vertex:(1 2 5 3 6 4 7 8),twist :(0 0 1 1 0 1 0 0)
R*U*R^(-1)*U*R*U^2*R^(-1)
 #vertex:(1 2 6 5 4 3 7 8),twist :(0 0 1 1 0 1 0 0)
R*F*U^2*F^(-1)*R^2*F*R
 #vertex:(1 2 5 3 6 4 7 8),twist :(0 0 1 0 1 1 0 0)
R*F*U^(-1)*F*R*F^2*R^(-1)
 #vertex:(1 2 6 5 4 3 7 8),twist :(0 0 1 0 1 1 0 0)
R*U^(-1)*R^(-1)*F*R^(-1)*F^(-1)*R
 #vertex:(1 2 6 5 3 4 7 8),twist :(0 0 1 1 0 1 0 0)
R*F^(-1)*U^(-1)*F*R^(-1)*U^(-1)*R
 #vertex:(1 2 5 6 4 3 7 8),twist :(0 0 1 0 1 1 0 0)
F*R*F^2*R^(-1)*U^2*R*F
 #vertex:(1 2 5 3 6 4 7 8),twist :(0 0 0 1 1 1 0 0)
F*R*F^(-1)*U*F*R^2*F^(-1)
 #vertex:(1 2 6 5 4 3 7 8),twist :(0 0 0 1 1 1 0 0)
F*U*F^2*U^(-1)*F^2*U*F
 #vertex:(1 2 5 3 6 4 7 8),twist :(0 0 1 1 1 0 0 0)
F*U*F^(-1)*U*F*U^2*F^(-1)
 #vertex:(1 2 6 5 4 3 7 8),twist :(0 0 1 1 1 0 0 0)
F*R^(-1)*F^(-1)*R*F^(-1)*U^(-1)*F
 #vertex:(1 2 4 5 6 3 7 8),twist :(0 0 0 1 1 1 0 0)
F*U^(-1)*F^(-1)*R*U^(-1)*R^(-1)*F
 #vertex:(1 2 6 3 4 5 7 8),twist :(0 0 1 1 1 0 0 0)
R^(-1)*U^2*R*U*R^(-1)*U*R
 #vertex:(1 2 6 5 4 3 7 8),twist :(0 0 0 1 1 1 0 0)
R^(-1)*F^2*R*U*R^(-1)*F*R
 #vertex:(1 2 6 5 4 3 7 8),twist :(0 0 1 1 1 0 0 0)
F^(-1)*R^2*F*R*U^(-1)*R*F
 #vertex:(1 2 6 5 4 3 7 8),twist :(0 0 1 0 1 1 0 0)
F^(-1)*U^2*F*U*F^(-1)*U*F
 #vertex:(1 2 6 5 4 3 7 8),twist :(0 0 1 1 0 1 0 0)
NIL
```

訳あって、ここまでで投稿します。後で、説明とソースを追加します。


[^1]: 通常のルービップキューブの解法解説では$R'$と表記しています。
[^2]: [群論の味わい](https://www.kyoritsu-pub.co.jp/bookdetail/9784320019416)では、「移動先において基準参照印に合わせるためひねる角度」を考えてます。
[^3]: $R$を定義すれば、残りの二つは計算できそうな気はしますが。
[^4]: もっとスマートな方法はあるはずですが。
[^5]: 値を持たないハッシュ表というのは邪道か？