LoginSignup
1
3

More than 3 years have passed since last update.

@futakuchi0117(田林 伸一郎)

１０次元立方体を眺める

Last updated at Posted at 2020-04-05

pythonで１０次元立方体を２次元に射影して回転アニメーションを生成します．

結果

コード

（矢印を押すと展開します）
import numpy as np
import matplotlib.pyplot as plt
import matplotlib.animation as animation
import matplotlib.collections as mc
# from mpl_toolkits.mplot3d import Axes3D 
import itertools
import time

def generate_rotation_2d(dim=3, axis1=0, axis2=1):
    """回転行列を生成する

    Args:
        dim(int): 次元
        axis1(int): 回転面の軸のインデックス
        axis2(int): 回転面の軸のインデックス

    Returns:
        function: 回転行列
    """
    def rotation_2d(theta):
        """回転行列

        Args:
            theta(float): 回転角

        Returns:
            np.2darray: 回転行列
        """
        rotation_matrix = np.zeros((dim, dim))
        for i in range(dim):
            rotation_matrix[i,i] = 1
        rotation_matrix[axis1, axis1] = np.cos(theta)
        rotation_matrix[axis1, axis2] = -np.sin(theta)
        rotation_matrix[axis2, axis1] = np.sin(theta)
        rotation_matrix[axis2, axis2] = np.cos(theta)
        return rotation_matrix
    return rotation_2d


def prod_rotation(theta, rotation_2d_list, dim):
    """回転行列の積をとる"""
    # 単位行列を生成
    rotation_matrix = np.identity(dim)
    # すべての回転行列の積をとる
    for rotation in rotation_2d_list:
        rotation_matrix = rotation_matrix @ rotation(theta)
    return rotation_matrix

def plot_cube(dim=3):
    """立方体を図示する"""
    # 回転行列の全体
    rotation_2d_list = [generate_rotation_2d(dim, axis1, axis2) for axis1, axis2 in itertools.combinations(range(dim),2)]

    # 頂点全体を生成
    vertex_list = [np.array(list(coordinate)).reshape(-1,1) for coordinate in itertools.product(*[[-1/2,1/2] for _ in range(dim)])]

    # エッジのインデックス全体を生成
    edge_index_list = []
    for comb in itertools.combinations(range(len(vertex_list)),2):
        if np.sum(np.abs(np.array(vertex_list[comb[0]]) - np.array(vertex_list[comb[1]]))) == 1:
            edge_index_list.append(comb)

    # 速度（設定用）
    speed = 4/dim
    # 回転角
    theta = (np.pi /180) * speed 

    # 回転行列を生成
    rotation_matrix = prod_rotation(theta, rotation_2d_list, dim)

    # 回転角で何回回転したかカウント
    cnt = -1

    def update(frame):
        plt.cla()

        nonlocal vertex_list, edge_index_list, rotation_matrix, dim, theta, cnt

        if cnt > 0:
            vertex_list = [rotation_matrix @ v for v in vertex_list]

        x = np.concatenate([v[0] for v in vertex_list])
        y = np.concatenate([v[1] for v in vertex_list])
#         z = np.concatenate([v[2] for v in vertex_list_rotated])

        lim = 1. + 0.1 * dim

        ax.set_xlim([-lim, lim])
        ax.set_xlabel("X")
        ax.set_ylim([-lim, lim])
        ax.set_ylabel("Y")
        ax.set_axis_off()

        ax.scatter(x,y, s=40/dim, c="red", alpha=0.6)

        # 線のリスト．　[(x0, y0), (x1, y1)]が1つの線
        lines = [[vertex_list[i1][:2].flatten(), vertex_list[i2][:2].flatten()] for i1, i2 in edge_index_list]
        colors = ["blue" for i in range(len(edge_index_list))]

        plt.text( 
            0.96, 0.96, f"$\Theta$: {theta*cnt*(180/np.pi)%360.0:.2f}[deg]",
            ha='right', va='top',
            transform=ax.transAxes
        )

        # 複数の線を追加
        lc = mc.LineCollection(lines, colors=colors, linewidths=1.0/dim)
        ax.add_collection(lc)

        cnt += 1

    fig, ax = plt.subplots(figsize=(5,5))
    ani = animation.FuncAnimation(fig, update, interval = 20, frames = int((np.pi*4)/theta)+1)
    ani.save(f'cube/cube_{dim}d_to_2d_720.mp4', writer="ffmpeg", dpi=145) 

if __name__ == "__main__":
    for dim in range(2,11):
        start = time.time()
        plot_cube(dim=dim)
        elapsed_time = time.time() - start
        print(f"dim: {dim}, elapsed_time:{elapsed_time}[sec]")

数学的なうんちく

２次元の回転行列は

\begin{pmatrix}
\cos\theta & -\sin\theta \\
\sin\theta & \cos\theta
\end{pmatrix}

です．

3次元の回転は

\begin{pmatrix}
\cos\theta & -\sin\theta & 0 \\
\sin\theta & \cos\theta & 0 \\
0 & 0 & 1
\end{pmatrix} , \qquad 

\begin{pmatrix}
\cos\theta & 0 & -\sin\theta \\
0 & 1 & 0 \\
\sin\theta & 0 & \cos\theta 
\end{pmatrix} , \qquad 
\begin{pmatrix}
1 & 0 & 0 \\
0 & \cos\theta & -\sin\theta \\
0 & \sin\theta & \cos\theta
\end{pmatrix}

の組み合わせで表現できます．

なので，自然な拡張として，N次元の回転を以下の行列の組み合わせで表現することができます：

\begin{align}
& \qquad \qquad \small{\text{i-th}} \qquad \qquad \small{\text{j-th}} \\
& \begin{pmatrix}
1 & \cdots & 0 & \cdots & 0 & \cdots &  0 \\
\vdots & \ddots & \vdots & & \vdots && \vdots \\
 0 & \cdots & \cos\theta & \cdots & -\sin\theta & \cdots & 0 \\
\vdots & & \vdots & & \vdots && \vdots \\
 0 & \cdots & \sin\theta & \cdots & \cos\theta & \cdots & 0 \\
\vdots & & \vdots & & \vdots &\ddots & \vdots \\
0 & \cdots & 0 & \cdots & 0 & \cdots &  1
\end{pmatrix} \begin{matrix} \\ \\  \small{\text{i-th}} \\ \\ \small{\text{j-th}}  \\ \\ \\ \end{matrix} ,  \qquad 1\le i < j \le N .
\end{align}

これを実装したのが上のコードです．

課題

１４次元までは上のコードで30分程度の計算時間で済んだのですが，それ以上次元を上げようとすると処理が終わらなくなりました．

処理が重くなりすぎたのは，回転行列による積の演算ではなく，頂点と辺のプロットでした．

N次元立方体の頂点の数は$2^N$，辺の数は$2^{N-1}N$です．

頂点の数は，10次元で1024, 20次元で1048576, 30次元で1073741824.

辺の数は，10次元で5120, 20次元で10485760, 30次元で16106127360.

恐ろしい増えっぷりですね．

その他の次元（観賞用）

無限次元へ

本記事で行ったプロットから察するに，無限次元立方体も２次元に射影すると，ほとんどいたる角度から見たとき「だいたい円盤」っぽく見えるんだろうとは思うのですが，実際のところはどうなんでしょうね？

高次元は不思議なことが多いのですが，それを遥かに超えた無限次元はまったく想像がつきません．

1
3
0

Register as a new user and use Qiita more conveniently

  1. You get articles that match your needs
  2. You can efficiently read back useful information
  3. You can use dark theme
What you can do with signing up
Sign upLogin
1
3