Eden モデルのアニメーション

境界面の成長モデルの可視化

---

以前の投稿でEden モデルという境界面の拡張等の研究で用いられているモデルを紹介した．

今回はそのモデルを理解しやすくしたいと思い，アニメーションを作成したので公開する．

コードは下記である．

import numpy as np
import matplotlib.pyplot as plt
from matplotlib.animation import FuncAnimation, PillowWriter
import io

def culture_spreading_model(size, steps, source_x, source_y):
    grid = np.zeros((size, size), dtype=int)
    age = np.zeros((size, size), dtype=int)
    grid[source_y, source_x] = 1
    current_culture = 1
    
    history = []
    
    for step in range(steps):
        occupied = np.where(grid > 0)
        neighbors = set()
        for i, j in zip(occupied[0], occupied[1]):
            neighbors.update([
                (i-1, j), (i+1, j),
                (i, j-1), (i, j+1)
            ])
        neighbors = [(i, j) for i, j in neighbors if 0 <= i < size and 0 <= j < size and grid[i, j] == 0]
        
        if neighbors:
            i, j = neighbors[np.random.randint(len(neighbors))]
            grid[i, j] = current_culture
            age[i, j] = step
        
        if np.random.random() < 0.03:
            current_culture += 1
            grid[source_y, source_x] = current_culture
            age[source_y, source_x] = step
        
        if step % 100 == 0:  # 100ステップごとに状態を保存
            history.append(grid.copy())
    
    return grid, age, history

# シミュレーションのパラメータ
grid_size = 100  # サイズを小さくして処理を軽くする
growth_steps = 10000
source_x, source_y = grid_size // 2, grid_size // 2

# シミュレーションの実行
result, age_map, history = culture_spreading_model(grid_size, growth_steps, source_x, source_y)

# アニメーションの作成
fig, ax = plt.subplots(figsize=(8, 8))

def update(frame):
    ax.clear()
    im = ax.imshow(history[frame], cmap='viridis', animated=True)
    ax.plot(source_x, source_y, 'ro', markersize=5)
    ax.set_title(f'Step: {frame * 100}')
    return [im]

ani = FuncAnimation(fig, update, frames=len(history), interval=100, blit=True)

# GIFとして保存
writer = PillowWriter(fps=10)
ani.save("culture_spreading.gif", writer=writer)

plt.close(fig)

# 結果の可視化（静的な図）
fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(16, 8))

im1 = ax1.imshow(result, cmap='viridis')
ax1.set_title('Final Culture Distribution')
fig.colorbar(im1, ax=ax1, label='Culture ID')
ax1.plot(source_x, source_y, 'ro', markersize=5)

im2 = ax2.imshow(age_map, cmap='YlOrRd')
ax2.set_title('Final Culture Age Distribution')
fig.colorbar(im2, ax=ax2, label='Age')
ax2.plot(source_x, source_y, 'bo', markersize=5)

plt.tight_layout()

# プロットをバイト列として保存
buffer = io.BytesIO()
plt.savefig(buffer, format='png')
buffer.seek(0)
plt.close(fig)

# 同心円状の分布を分析
def analyze_concentric_distribution(grid, center_x, center_y):
    max_radius = min(grid_size - center_x, center_x, grid_size - center_y, center_y)
    culture_counts = {}
    
    for r in range(1, max_radius):
        cultures = set()
        for theta in np.linspace(0, 2*np.pi, 100):
            x = int(center_x + r * np.cos(theta))
            y = int(center_y + r * np.sin(theta))
            cultures.add(grid[y, x])
        culture_counts[r] = len(cultures)
    
    return culture_counts

concentric_distribution = analyze_concentric_distribution(result, source_x, source_y)

# 同心円状の分布をプロット
plt.figure(figsize=(10, 5))
plt.plot(list(concentric_distribution.keys()), list(concentric_distribution.values()))
plt.title('Number of Distinct Cultures vs Distance from Center')
plt.xlabel('Distance from Center')
plt.ylabel('Number of Distinct Cultures')

# プロットをバイト列として保存
buffer2 = io.BytesIO()
plt.savefig(buffer2, format='png')
buffer2.seek(0)
plt.close()

print("シミュレーションが完了し、結果が保存されました。")
print("1. 'culture_spreading.gif': 文化伝播のアニメーション")
print("2. 最終的な文化分布と年齢分布の図が生成されました。")
print("3. 同心円状の分布分析の図が生成されました。")

実行結果は下記．