はじめに
SimPyというPythonの離散事象シミュレーション用のパッケージを見つけて試してみたら気に入ったので自分用の備忘録も兼ねて使い方をまとめていく．
- SimPyのドキュメント
- SimPyのソース
- 第1回の記事（はじめの一歩）
- 第2回の記事（リソースを理解しよう：Container編）
- 第3回の記事（リソースを理解しよう：StoreとResource編）
第4回目の今回は，シミュレーションの途中経過のアニメーションに挑戦しよう．といっても，SimPy自体にはアニメーションの機能は存在しないので別の手段の助けを借りることになる．Pythonにはさまざまなビジュアライゼーションのツールが存在するので，基本的には好みのものを利用すればよい．ここでは，それらのうち，
の2つの手段を用いた例をそれぞれ簡単に紹介することにする．
matplotlibを用いた例
第2回に取り上げた在庫管理モデルの例を思い出そう．Containerリソースをデフォルトのまま用いた例である．リソースの属性levelを参照することで在庫量がわかるので，その推移を折れ線グラフで表示してみる．
import numpy as np
import matplotlib.pyplot as plt
import random
import simpy
def manager(env):
env.model.ordered = False # no back order to receive
env.stocktake = env.event() # create the first signal (event)
while True:
yield env.stocktake
report(env)
if not env.model.ordered and env.model.level <= 10:
# only when no back order to receive
# reorder point = 10
env.process(deliverer(env)) # activate deliverer
env.model.ordered = True # back order will be received
env.stocktake = env.event() # create the next signal (event)
def deliverer(env):
yield env.timeout(5) # delivery lead time = 5
env.model.put(20) # back order is recieved
env.model.ordered = False # no back order to receive
def customer(env):
while True:
time_to = random.expovariate(1)
yield env.timeout(time_to)
how_many = random.randint(1, 3)
env.model.get(how_many)
env.stocktake.succeed() # signal for stocktaking (event)
def report(env):
print('[{}] current level: {}, orderd: {}, queue length: {} '.format(round(env.now), env.model.level, env.model.ordered, len(env.model.get_queue)))
# ---------- code for visualization ----------
env.x.append(env.now) # time
env.y.append(env.model.level) # container's level
# ---------- ---------- ---------- ----------
def main():
env = simpy.Environment()
env.model = simpy.Container(env, init=10) # model is marely a Container
env.process(manager(env))
env.process(customer(env))
# env.run(until=200)
# ---------- code for visualization ----------
env.x = [0] # time
env.y = [10] # container's level
fig = plt.figure(1)
ax = fig.add_subplot(111)
ax.set_xlabel('time')
ax.set_ylabel('number')
ax.set_xlim(0, 200)
ax.set_ylim(0, 30)
line, = ax.plot(env.x, env.y, label='level')
ax.legend()
ax.grid()
for t in range(1, 200):
env.run(until=t) # stepwise execution
line.set_data(env.x, env.y)
plt.pause(0.1)
plt.show()
# ---------- ---------- ---------- ----------
if __name__ == "__main__":
main()
main()関数の後半に折れ線グラフを表示するコードが追加されている．ポイントは，for文で時刻tを1から200まで1ずつ増加させ，各ループでenv.run(until=t)を呼んでいることである．これによって時間1刻みでシミュレーションを進め，そのたびに折れ線グラフを更新していっている．
なお，表示に使うデータ（env.x，env.y）は，このfor文の中で更新してもいいが，時間1よりも細かい単位での変化を反映させるために，上の例では，関数report()の中で更新している．
p5（p5py）を用いた例
上の例のように，折れ線グラフなどのグラフ表示でシミュレーションの経緯を追いかけたい場合は，matplotlibが便利である．しかし，ときには，もう少し自由度の高いアニメーションを作成したくなることもある．その場合は，いわゆるドロー系のツールに頼りたい．
ここでは，Pythonで利用できるドロー系のパッケージの一例として，p5（p5py）を利用してみる．これは，ProcessingをPythonに移植したものだと考えればいいだろう．もしかするとあまりメジャーなパッケージではないかもしれないが，（筆者のような）Processingユーザにとってはありがたい選択肢である．
p5（p5py）とProcessingの違いはここにまとめられている．なお，余談として，Processing自体にもPython Modeとよばれる機能があるが，これはJythonなので，SimPyと組み合わせるのは難しい．
p5（p5py）の基本的な構文は下のようになる．
from p5 import *
def setup():
size(**スケッチの幅**, **スケッチの高さ**)
# setup()関数の中に初期設定コードを書く
def draw():
background(**背景色の指定**)
# draw()関数の中に描画のコードを書く
if __name__ == '__main__':
run() # 実際に描画する際にrun（）関数を呼ぶ
draw()関数はループ実行されるので，これをアニメーションのフレームに対応させる．最後に，描画のためのp5pyのrun()関数を呼ぶ必要があるので，SimPyの方のrun()関数はdraw()関数の中でアニメーションのペースに合わせて少しずつ実行していくようにするといいだろう．
例をみてみよう．
import random
import simpy
from p5 import *
def manager(env):
env.model.ordered = False # no back order to receive
env.stocktake = env.event() # create the first signal (event)
while True:
yield env.stocktake
report(env)
if not env.model.ordered and env.model.level <= 10:
# only when no back order to receive
# reorder point = 10
env.process(deliverer(env)) # activate deliverer
env.model.ordered = True # back order will be received
env.stocktake = env.event() # create the next signal (event)
def deliverer(env):
yield env.timeout(5) # delivery lead time = 5
env.model.put(20) # back order is recieved
env.model.ordered = False # no back order to receive
def customer(env):
while True:
time_to = random.expovariate(1)
yield env.timeout(time_to)
how_many = random.randint(1, 3)
env.model.get(how_many)
env.stocktake.succeed() # signal for stocktaking (event)
def report(env):
print('[{}] current level: {}, orderd: {}, queue length: {} '.format(round(env.now), env.model.level, env.model.ordered, len(env.model.get_queue)))
def simpy_setup(): # setup of simpy environment
env = simpy.Environment()
env.model = simpy.Container(env, init=10) # model is marely a Container
env.process(manager(env))
env.process(customer(env))
# env.run(until=200)
return env
# ---------- code for visualization ----------
def setup(): # setup of p5py sketch
size(600, 600)
def draw(): # drawing loop of visualization by p5py
background(255)
if frame_count <= 200:
simpy_env.run(until=frame_count) # keep pace with draw()
level = simpy_env.model.level
queue = len(simpy_env.model.get_queue)
ordered = simpy_env.model.ordered
else:
exit()
translate(250, 350)
line((-200, 0), (300, 0))
if queue > 0:
fill(127, 0, 0)
rect((0,0), 100, 10 *queue) # how many customers are waiting
if level > 0:
fill(127)
rect((0,0), 100, -10 *level) # how many items at hand
translate(0, -10 *level)
if ordered:
fill(255)
rect((0,0), 100, -10 *20) # size of back order to arrive
# ---------- ---------- ---------- ----------
if __name__ == "__main__":
simpy_env = simpy_setup()
run(frame_rate=10) # you can control the speed here
もとの例にあったmain()関数を微修正した，simpy_setup()という名称の関数を用意した．これは，SimPy側の設定を担当する関数であり，環境envを返す．この戻り値をsimpy_envという変数で受け取って，p5py側で利用している．
アニメーションの内容は実際に実行して確認してみてもらえたらと思うが，在庫量とバックオーダーが積み上げ型の棒グラフで表示され，欠品が出たときは，待機中の顧客の人数が下方向の棒グラフでわかるようになっている．
まとめ
今回はSimPyで実行しているシミュレーションの途中経過をアニメーション表示する方法の例を紹介した．ここで紹介した以外にも様々な方法が考えられると思う．コンソールに文字が表示されるだけでは味気ないので，ぜひアニメーションにも挑戦してみてほしい．次回はシミュレーションにエージェントを追加する例を簡単に紹介したい．