ガウゼの法則のシミュレーション

ガウゼの法則をPythonで学ぶ

---

以前の投稿でガウゼの法則について見た．

今回はこのガウゼの法則をシミュレーションで可視化した．

コードは下記．

import numpy as np
import matplotlib.pyplot as plt
import matplotlib.colors as mcolors
import imageio

# パラメータの設定
GRID_SIZE = 50  # グリッドのサイズ
num_agents_A = 500  # エージェントAの初期数
num_agents_B = 500  # エージェントBの初期数
p_reproduce_A = 0.05  # エージェントAの繁殖確率
p_reproduce_B = 0.03  # エージェントBの繁殖確率（Aよりわずかに低い）
p_die = 0.02  # 死亡確率
num_time_steps = 300  # シミュレーションのステップ数

# グリッドの初期化
grid = np.zeros((GRID_SIZE, GRID_SIZE), dtype=int)
positions = np.random.choice(GRID_SIZE * GRID_SIZE, num_agents_A + num_agents_B, replace=False)
positions_A = positions[:num_agents_A]
positions_B = positions[num_agents_A:]

grid.flat[positions_A] = 1  # エージェントAを配置
grid.flat[positions_B] = 2  # エージェントBを配置

# 隣接セルを取得する関数
def get_neighbors(i, j, grid_size):
    neighbors = []
    for x in [-1, 0, 1]:
        for y in [-1, 0, 1]:
            if x == 0 and y == 0:
                continue
            ni, nj = i + x, j + y
            if 0 <= ni < grid_size and 0 <= nj < grid_size:
                neighbors.append((ni, nj))
    return neighbors

# カラーマップの設定（0: 白、1: 青、2: 赤）
cmap = mcolors.ListedColormap(['white', 'blue', 'red'])

# シミュレーションとGIFの作成
with imageio.get_writer('gause_simulation.gif', mode='I', fps=10) as writer:
    for t in range(num_time_steps):
        new_grid = grid.copy()

        for i in range(GRID_SIZE):
            for j in range(GRID_SIZE):
                if grid[i, j] == 0:
                    continue
                agent_type = grid[i, j]
                # 死亡判定
                if np.random.rand() < p_die:
                    new_grid[i, j] = 0
                    continue
                # 繁殖判定
                p_reproduce = p_reproduce_A if agent_type == 1 else p_reproduce_B
                if np.random.rand() < p_reproduce:
                    neighbors = get_neighbors(i, j, GRID_SIZE)
                    np.random.shuffle(neighbors)
                    for ni, nj in neighbors:
                        if new_grid[ni, nj] == 0:
                            new_grid[ni, nj] = agent_type
                            break
                # 移動判定
                if np.random.rand() < 0.5:
                    neighbors = get_neighbors(i, j, GRID_SIZE)
                    np.random.shuffle(neighbors)
                    for ni, nj in neighbors:
                        if new_grid[ni, nj] == 0:
                            new_grid[ni, nj] = agent_type
                            new_grid[i, j] = 0
                            break
        grid = new_grid.copy()
        # グリッドのプロット
        fig, ax = plt.subplots(figsize=(6, 6))
        ax.imshow(grid, cmap=cmap, vmin=0, vmax=2)
        ax.axis('off')
        plt.tight_layout()
        # フレームをGIFに追加
        fig.canvas.draw()
        image = np.frombuffer(fig.canvas.tostring_rgb(), dtype='uint8')
        image = image.reshape(fig.canvas.get_width_height()[::-1] + (3,))
        writer.append_data(image)
        plt.close()