やること
ライフゲームは1ステップ先を計算するのは簡単ですが、1ステップ前を計算するのはなかなか難しいです。vcoptで逆ライフゲームを試してみましょう。今回はvcopt().dcGA()を繰り返し実行するのでトリッキーです。
参考文献
ライフゲームは、順方向の結果は一意に定まりますが、逆方向は一意ではありません。したがって、ハッシュや暗号に応用できると言われています。

https://www.jstage.jst.go.jp/article/pjsai/JSAI08/0/JSAI08_0_196/_pdf
vcoptの使い方についてはチュートリアルをご参照ください。
vcoptの仕様については最新の仕様書をご参照ください。本記事執筆時とは仕様が異なる場合があります。
実行環境
WinPython3.6をおすすめしています。

WinPython - Browse /WinPython_3.6/3.6.7.0 at SourceForge.net
Portable Scientific Python 2/3 32/64bit Distribution for Windows
ライフゲームの基本
こちらの記事をご参照ください。
ソースコード
長い。
import numpy as np
from scipy import signal
import matplotlib.pyplot as plt
from vcopt import vcopt
#============================
#設定
#============================
#高さ、幅
h, w = 10, 10
#任意の最終状態を用意
state = np.array([[0,0,1],
[1,0,1],
[0,1,1]])
#フィールドのどこに置くか(左上点を指定)
p = (4, 4)
#============================
#関数など
#============================
#目標フィールドの生成
f_aim = np.zeros((h, w), dtype=bool)
#任意の状態を置く
f_aim[p[0]:p[0]+len(state), p[1]:p[1]+len(state[0])] = state
#畳み込み用のフィルタ
g = np.array([[1.0, 1.0, 1.0],
[1.0, 0.0, 1.0],
[1.0, 1.0, 1.0]])
h_crop, w_crop = None, None
p_crop = (None, None)
#paraからfを復元する関数
def para_to_f(para, h_crop, w_crop, p_crop):
#クロップ部分を復元
f_crop = para.reshape(h_crop, w_crop)
#フィールド全体の復元
f = np.zeros((h, w), dtype=bool)
f[p_crop[0]:p_crop[0]+h_crop, p_crop[1]:p_crop[1]+w_crop] = f_crop
return f
#fからfutureを計算する関数
def f_to_future(f):
#周囲の生存マスを足し込む
mask = signal.convolve2d(f, g, mode='same', boundary='wrap')
#1ステップ未来のフィールド(すべて死状態)
future = np.zeros(f.shape, dtype=int)
#生きているマスが生きる条件(=生存)
future[mask*f==2] = 1
future[mask*f==3] = 1
#死んでいるマスが生きる条件(=誕生)
future[mask*~f==3] = 1
return future
#3つのグラフを並べて表示する関数
def show_3(f1, f2, f3, title1, title2, title3, i):
plt.figure(figsize=(10, 10))
plt.subplot(1, 3, 1)
plt.imshow(f1, cmap='inferno')
plt.title(title1)
plt.subplot(1, 3, 2)
plt.imshow(f2, cmap='inferno')
plt.title(title2)
plt.subplot(1, 3, 3)
plt.imshow(f3, cmap='inferno')
plt.title(title3)
plt.savefig('save/{}.png'.format(i), bbox_inches='tight', pad_inches=0)
plt.show(), print()
#生きているセルのある部分を余白つけてクロップする関数
def crop(f):
col_sum = np.sum(f, axis=0)
row_sum = np.sum(f, axis=1)
col_exist = np.where(col_sum>0)[0]
col_start = max(col_exist[0] - 1, 0)
col_end = min(col_exist[-1] + 2, f.shape[1])
row_exist = np.where(row_sum>0)[0]
row_start = max(row_exist[0] - 1, 0)
row_end = min(row_exist[-1] + 2, f.shape[0])
f_crop = f[row_start:row_end, col_start:col_end]
h_crop, w_crop = f_crop.shape[0], f_crop.shape[1]
p_crop = (row_start, col_start)
return f_crop, h_crop, w_crop, p_crop
#評価関数
def rGL_score(para):
#paraからfの復元
f = para_to_f(para, h_crop, w_crop, p_crop)
#1ステップ未来のフィールド
future = f_to_future(f)
#ステップ前後の差分マス数
diff = np.sum((future - f_aim)**2)
#前ステップの生セル比率
l_score = np.sum(f) / (h * w)
#スコアを返す
return diff + l_score
#============================
#リバース開始
#============================
score = 0.0
step = -1
while score < 1.0:
print('{}stepを計算中'.format(step))
#生きているセルのある部分を余白つけてクロップ
f_crop, h_crop, w_crop, p_crop = crop(f_aim)
#パラメータ範囲
para_range = [[i for i in range(2)] for j in range(h_crop * w_crop)]
#GAで最適化
para, score = vcopt().dcGA(para_range,
rGL_score,
0.0,
show_pool_func='bar',
seed=1,
pool_num=1000)
#結果の表示
print('score={}'.format(score))
#paraからfの復元
f = para_to_f(para, h_crop, w_crop, p_crop)
#1ステップ未来のフィールド
future = f_to_future(f)
#表示
show_3(f, future, f_aim, 'before', 'after', 'f_aim', -step)
#更新
f_aim = f
step -= 1
評価関数を少し工夫しています。あるparaからフィールドを復元し、1ステップ後を計算し、目標フィールドとの差分をdiffとします。これだけだとセルが拡散しがちなので、セルがまとまりやすいようにゴニョゴニョしています(説明が面倒くさくなり絶命)。
結果
ソースコードをそのまま実行すると、グライダーの逆ライフゲーム探索が始まります。
-1stepを計算中
[|< ] gen=73000, best_score=0.05
score=0.05


-2stepを計算中
[ |< ] gen=55000, best_score=0.05
score=0.05


-3stepを計算中
[|< ] gen=87000, best_score=0.05
score=0.05


-4stepを計算中
[|< ] gen=77000, best_score=0.06
score=0.06


-5stepを計算中
[|< ] gen=71000, best_score=0.08
score=0.08


-6stepを計算中
[|< ] gen=143000, best_score=0.1
score=0.1


-7stepを計算中
[|< ] gen=138000, best_score=0.33
score=0.33


-8stepを計算中
[ |< ] gen=194000, best_score=3.23
score=3.23


ということで、-7ステップ目で急激に発散し、-8ステップ目は見つかりませんでした。score=3.23というのは、3マスの食い違いが生じ、フィールドの23%が生セルという意味です。
順方向で確認
-7ステップ目から再生してみます。
import numpy as np
from scipy import signal
import matplotlib.pyplot as plt
#============================
#初期状態の設定
#============================
#高さ、幅
h, w = 10, 10
#任意の状態を用意
state = np.array([[0,0,1,0,1,1,1],
[0,1,1,0,1,0,0],
[1,0,1,1,0,0,0],
[0,0,1,1,1,1,1],
[0,1,0,0,0,1,0],
[0,0,1,1,0,1,0],
[0,1,1,0,1,0,0],
[0,0,1,0,1,1,0],
[1,0,1,1,0,1,0],
[0,0,1,1,1,0,0]])
#フィールドのどこに置くか(左上点を指定)
p = (0, 1)
#終了ステップ数
max_step = 15
#============================
#メイン処理
#============================
#フィールドの生成
f = np.zeros((h, w), dtype=bool)
#任意の状態を置く
f[p[0]:p[0]+len(state), p[1]:p[1]+len(state[0])] = state
#畳み込み用のフィルタ
g = np.array([[1.0, 1.0, 1.0],
[1.0, 0.0, 1.0],
[1.0, 1.0, 1.0]])
#初期状態の表示
#plt.figure(figsize=(10, 10))
plt.imshow(f, cmap='inferno')
#plt.savefig('save/{}.png'.format(0), bbox_inches='tight', pad_inches=0)
plt.show(), print()
#状態の更新
for i in range(1, max_step + 1):
#周囲の生存マスを足し込む
mask = signal.convolve2d(f, g, mode='same', boundary='wrap')
#1ステップ未来のフィールド(すべて死状態)
future = np.zeros(f.shape, dtype=bool)
#生きているマスが生きる条件(=生存)
future[mask*f==2] = 1
future[mask*f==3] = 1
#死んでいるマスが生きる条件(=誕生)
future[mask*~f==3] = 1
#フィールドの更新(浅いコピーに注意)
f = future
#表示
#plt.figure(figsize=(10, 10))
plt.imshow(f, cmap='inferno')
#plt.savefig('save/{}.png'.format(i), bbox_inches='tight', pad_inches=0)
plt.show(), print()
結果

