pythonのsimpyでシミュレーションをしてみた その３

今回の記事は[前回](http://qiita.com/ogata-k/items/efa963d5752a64490de2)の予告通りの総集編です。この記事でこのシリーズは最後にしようと思います。
#今回するシミュレーション
今回解説するシミュレーションは一言でいうと商品が選択できるかのシミュレーションです。詳しく説明するとほしい商品をほしい個数だけ買う客がいるとします。この客は買い物をするときに商品棚に並んで商品を選んでいくと思います。また店側からすれば商品を補充して商品を買わない客が出ないようにしたいと思います。その割合を調べるためのシミュレーションです。
#コードの全体像
最初にコードの全体を載せます。少し例としては長いかもしれません。あとで解説をしますがコメントを見ればなんとなくはわかるかもしれません。

```python:ShelfStore.py
# -*- coding: utf-8 -*-

import simpy
import random
from scipy.stats import poisson
import matplotlib.pyplot as plt
import sys
import os


RAND_ITEM = [1, 4]  # [min, max]顧客の商品選択数の乱数用
CUSTOMER_MU = 5  # 顧客製造のpoisson分布用の定数

SHELF_MAX = 15  # 棚に並べる商品数の最大値
SHELF_INIT = 5  # 棚に並べてある商品数の初期値
SHELF_INTERVAL = 15  # 棚に並べる時間の間隔
SHELF_LIMIT_NUM = 5  # 商品数の下限値
SHELF_MU = 18  # 商品の補充タイミングのpoisson分布用の定数
SHELF_NUM_ITEMS = 4  # 商品を並べる個数の一定値
SHELF_MAX_CUSTOMER = 3  # 棚に並ぶことのできる人数の限界
PUT_SPEED = 1  # 商品一つあたりに並べるのにかかる時間

choice = 0  # 選ぶことのできた人の数
choices = []  # たくさんのシミュレーションの結果のchoiceを収容するリスト
no_choice = 0  # 選ぶことのできなかった人の数
no_choices = []  # たくさんのシミュレーションの結果のno_choiceを収容するリスト

SIM_TIME = 100000  # シミュレートする時間

WHEN_PUT_FUNC_LIST = ["Interval", "Fewer", "Poisson"]  #補充タイミングの関数名の指定用
PUT_FUNC_LIST = ["Maximum", "Const"]  # 補充方法の関数名の指定用
texts = []  # シミュレートの実行時の条件


class Shelf:
    """
    商品棚のクラスで店側が商品を補充するようなタイミングを自由に設定できる。
    """
    def __init__(self, env, when_put_type="Interval", put_type="Maximum"):
        """
        when_put_typeは商品補充のタイミング用の関数指定\n
        put_typeは商品補充時の商品補充方法の関数指定
        """
        self.shelf = simpy.Resource(env, capacity=SHELF_MAX_CUSTOMER)
        self.items = simpy.Container(env, capacity=SHELF_MAX, init=SHELF_INIT)
        self.flag = 1  # 商品補充中でないなら1、補充中なら０となる
        env.process(self.Switch_Func(env, when_put_type, put_type))

    def Switch_Func(self, env, when_put_type="Interval", put_type="Maximum"):
        """シミュレート用の関数の切り替え用の関数"""
        if when_put_type == "Interval":
            when_put_func = self.PutInterval
        elif when_put_type == "Fewer":
            when_put_func = self.PutWhenFewer
        elif when_put_type == "Poisson":
            when_put_func = self.PutPoisson

        if put_type == "Maximum":
            put_func = self.PutMax
        elif put_type == "Const":
            put_func = self.PutConst

        return self.ChoiceWhenPutFunc(env, when_put_func, put_func)

    def ChoiceWhenPutFunc(self, env, when_put_func, put_func):
        """利用する関数を呼び出すための関数"""
        return when_put_func(env, put_func)

    def PutInterval(self, env, put_func):
        """決まった間隔で商品を補充する関数"""
        while True:
            yield env.timeout(SHELF_INTERVAL)
            if self.items.level < self.items.capacity and self.flag:
                self.flag = 0
                env.process(put_func(env))

    def PutWhenFewer(self, env, put_func):
        """少なくなったら商品を補充する関数"""
        print("%d個より少なくなったら商品を追加します。" % SHELF_LIMIT_NUM)
        while True:
            yield env.timeout(1)
            if self.items.level < SHELF_LIMIT_NUM and self.flag:
                self.flag = 0
                env.process(put_func(env))

    def PutPoisson(self, env, put_func):
        """poisson分布に従ってランダムにタイミングを決めて商品を補充する関数"""
        print("mu=%dのpoisson分布に従うようなタイミングで商品を追加します。" % SHELF_MU)
        while True:
            stats = poisson.rvs(SHELF_MU, size=1)
            yield env.timeout(stats)
            if self.items.level < self.items.capacity and self.flag:
                self.flag = 0
                env.process(put_func(env))

    def PutMax(self, env):
        """商品の補充時に商品棚の最大商品量になるように補充するような関数"""
        amount = self.items.capacity - self.items.level
        print("%.1fに%d個の商品を追加し始めました。"
              % (env.now, amount))
        yield env.timeout(amount*PUT_SPEED)
        yield self.items.put(amount)
        print("%.1fに商品の追加を終えて残りは%d個です。"
              % (env.now, self.items.level))
        self.flag = 1

    def PutConst(self, env):
        """
        決められた定数の個数だけ商品を補充する関数。
        もし商品補充後の商品数が商品棚の最大商品量を超えたら最大商品量になるように商品を補充
        """
        amount = SHELF_NUM_ITEMS
        if amount + self.items.level > self.items.capacity:
            amount = self.items.capacity - self.items.level
        print("%.1fに%d個の商品を追加し始めました。"
              % (env.now, amount))
        yield env.timeout(amount*PUT_SPEED)
        yield self.items.put(amount)
        print("%.1fに商品の追加を終えて残りは%d個です。"
              % (env.now, self.items.level))
        self.flag = 1


class Customer:
    """
    商品を選ぶ客を表現するクラス。商品の選ぶ個数はランダムだがほしい個数分の商品がなければあきらめて帰る。
    """
    def __init__(self, env, shelf):
        self.item = random.randint(*RAND_ITEM)  # 選択する商品の個数
        env.process(self.select(env, shelf))

    def select(self, env, shelf):
        """商品を選ぶ関数で商品を選べた個数と選べなかった個数を記録する関数。"""
        global choice, no_choice
        with shelf.shelf.request() as req:
            yield req
            if shelf.items.level < self.item:
                no_choice += 1
                print("%.1fに商品を選びに来たのですが%d個の商品がなくて帰ってしまいました。"
                     % (env.now, self.item))
            else:
                choice += 1
                yield shelf.items.get(self.item)
                print("%.1fに%d個の商品を選んでいきました。" % (env.now, self.item))
                print("商品棚に残っている商品は%d個です。" % shelf.items.level)


def customer_fuctory(env, shelf):
    """名前の通り客を作る工場用の関数。これによって疑似的な客の流れが表せる。"""
    while True:
        time = poisson.rvs(CUSTOMER_MU, size=1)
        yield env.timeout(time)
        Customer(env, shelf)


def percentage(a, b):
    """入力された数値の割合を返す関数。"""
    sigma = a+b
    return a/sigma, b/sigma


def main(when_put_func_number, put_func_number):
    """メイン関数"""
    global choice, no_choice
    choice = 0
    no_choice = 0
    env = simpy.Environment()
    shelf = Shelf(env, WHEN_PUT_FUNC_LIST[when_put_func_number], PUT_FUNC_LIST[put_func_number])
    print("商品補充タイミング：%s、商品補充方法：%s"
          % (WHEN_PUT_FUNC_LIST[when_put_func_number], PUT_FUNC_LIST[put_func_number]))
    env.process(customer_fuctory(env, shelf))
    env.run(until=SIM_TIME)
    print("\n結果\nchoice\tno_choice\t合計")
    sigma = choice+no_choice
    print("%.3f\t%.3f\t1.00" % percentage(choice, no_choice))
    print("%d人\t%d人\t%d人" % (choice, no_choice, sigma))
    print("商品補充タイミング：%s  商品補充方法：%s"
          % (WHEN_PUT_FUNC_LIST[when_put_func_number], PUT_FUNC_LIST[put_func_number]))

    return choice, no_choice, sigma

# resultフォルダをなければ作成してそのフォルダに移動
if not os.path.isdir("./result"):
    os.mkdir("./result")
os.chdir("./result")

"""ここのループではそれぞれの補充タイミングと補充方法を指定してそれぞれのファイルに出力"""
for i in range(len(WHEN_PUT_FUNC_LIST)):
    for j in range(len(PUT_FUNC_LIST)):
        f_name = "%sAnd%s.txt" % (WHEN_PUT_FUNC_LIST[i], PUT_FUNC_LIST[j])
        with open(f_name, "w") as fd:
            sys.stdout = fd
            c, nc, sigma = main(i, j)
        choices.append(c/sigma)
        no_choices.append(nc/sigma)
        texts.append(WHEN_PUT_FUNC_LIST[i] + " and " + PUT_FUNC_LIST[j] +
                     "\nnum=%d" % sigma)

n = len(WHEN_PUT_FUNC_LIST) * len(PUT_FUNC_LIST)  # 組み合わせの総数
# テキストとして描画
f_name = "all_result.txt"
with open(f_name, "w") as fd:
    sys.stdout = fd
    for i in range(n):
        print("\n結果")
        print(texts[i]+"(全体)")
        print("選択できた人の割合", end=" ")
        print(choices[i])
        print("選択できなかった人の割合", end=" ")
        print(no_choices[i])

# 描画
plt.figure(figsize=(20, 15))
left = [i for i in range(n)]
plt.bar(left, choices, color="green", bottom=no_choices, align="center",
        label="choice")
plt.bar(left, no_choices, color="blue", align="center",
        tick_label=texts, label="no_choice")
plt.legend()
plt.ylim(0, 1)
plt.xlim(-0.5, n-0.5)
plt.title("simulation of Shelf's Items\nsim time = %s" % SIM_TIME)
plt.xlabel("way of simulation")
plt.ylabel("percentage of customer")

plt.savefig("Shelf_Item_Sim.jpg")
# plt.show()
```

このコードではたくさんの標準出力があるためその出力先をresultフォルダに変更しています。そのままコピペで実行をしても何も出力されないのはその為なのでミスではありません。またresultフォルダが既に存在する状態でこのプログラムを実行するとそこに出力されるので気を付けてください。
# コード説明
## Shelfクラス
インスタンス化したときに商品棚のスペース（Resourceクラス）や所有している商品（Containerクラス）の定義の後、引数としてシミュレーションする補充タイミングと補充方法を指定した条件でenv環境のプロセスに関数を追加しています。何も指定しない場合は等間隔で最大量になるように補充する設定になっています。ここで補充するタイミングを指定する関数のうち少なくなったら補充する関数のPutWhenFewerを使って解説をします。

```
    def PutWhenFewer(self, env, put_func):
        """少なくなったら商品を補充する関数"""
        print("%d個より少なくなったら商品を追加します。" % SHELF_LIMIT_NUM)
        while True:
            yield env.timeout(1)
            if self.items.level < SHELF_LIMIT_NUM and self.flag:
                self.flag = 0
                env.process(put_func(env))

```

わかりやすいように関数の定義の部分を抜き出してみました。この関数は無限ループを用いた無限イテレータです。そしてenv環境で時間が進むたびにself.items.levelで商品数を取得してその数が商品数下限値を意味する定数より小さいかつ商品補充中でないときにプロセスに追加する仕様です。
次は補充方法の関数の例としてPutConst関数を説明します。


```
    def PutConst(self, env):
        """
        決められた定数の個数だけ商品を補充する関数。
        もし商品補充後の商品数が商品棚の最大商品量を超えたら最大商品量になるように商品を補充
        """
        amount = SHELF_NUM_ITEMS
        if amount + self.items.level > self.items.capacity:
            amount = self.items.capacity - self.items.level
        print("%.1fに%d個の商品を追加し始めました。"
              % (env.now, amount))
        yield env.timeout(amount*PUT_SPEED)
        yield self.items.put(amount)
        print("%.1fに商品の追加を終えて残りは%d個です。"
              % (env.now, self.items.level))
        self.flag = 1
```
この関数は最初にamount変数に商品補充の固定量を代入しています。その後今ある商品にamount個の商品補充したときに最大商品数を超えないようにするためamount変数を調整します。そして補充する商品数に応じた時間をかけて補充していきます。ここで先に商品の補充する要求のyield self.items.put(amount)を呼び出すと補充してから補充する時間待つことになるので商品の追加が終わる前に客が補充した商品を選択することができるようになってしまいます。
## Customerクラス

```
class Customer:
    """
    商品を選ぶ客を表現するクラス。商品の選ぶ個数はランダムだがほしい個数分の商品がなければあきらめて帰る。
    """
    def __init__(self, env, shelf):
        self.item = random.randint(*RAND_ITEM)  # 選択する商品の個数
        env.process(self.select(env, shelf))

    def select(self, env, shelf):
        """商品を選ぶ関数で商品を選べた個数と選べなかった個数を記録する関数。"""
        global choice, no_choice
        with shelf.shelf.request() as req:
            yield req
            if shelf.items.level < self.item:
                no_choice += 1
                print("%.1fに商品を選びに来たのですが%d個の商品がなくて帰ってしまいました。"
                     % (env.now, self.item))
            else:
                choice += 1
                yield shelf.items.get(self.item)
                print("%.1fに%d個の商品を選んでいきました。" % (env.now, self.item))
                print("商品棚に残っている商品は%d個です。" % shelf.items.level)


```

ここでもわかりやすく定義の部分を抜き出してみました。Customerクラスのインスタンス化の時の挙動もShelfクラスのインスタンス化の時の挙動と似たようなものです。ですのでselect関数の仕様のみを説明したいと思います。[前回](http://qiita.com/ogata-k/items/efa963d5752a64490de2)で説明したwith文を使って商品棚を要求します。商品棚の定義の時にcapacityとして定義した値を商品棚の個数ではなく商品棚を利用できる客の個数として定義したのがここに効いてきます。これができるのはContainerクラスが商品棚（Resourceクラス）のcapacityの個数一つ一つと対応しているのではなくenv環境に対応しているからです。そして商品棚を使うことができた客は商品を選ぶのですが、ここですぐにyield shelf.items.get(selh.item)としてしまうと客は商品を手に入れるまで待ってしまいます。この関数ではそうならないように場合分けをしています。その場合分けにより商品がほしい個数だけ商品棚にある客は少しの時間も待たずに商品を選んでいきます。そしてこの関数を抜けるのでイテレータとしてのエラーを返し、この関数は役目を終えます。
# まとめ
これまででsimpyの簡単な全体像を想像できる程度には説明ができたかと思います。simpyは今までの説明のように離散型の待ち行列を簡単にシミュレーションするすべを与えてくれます。さらにsimpyを用いればそれぞれのクラスごとに行動をまとめ、感覚的に記述していくことができます。行動順序さえ与えてしまえば動いてくれるのでシミュレーション初心者にとっても使いやすいライブラリになると思います。（筆者も初心者の一人ですが。）simpyはpythonライブラリの中ではあまり知る人が少ない部類に入る上、有名なsympyと名前がかぶっていて情報が探しにくく、使える状況が限られています。しかし３回にわたって説明したように強力で面白いライブラリなので使う人が増えてくれるとありがたいです。そして情報を増やしていってください。（切実）
以上でこのシリーズを終わりにします。

＊筆者は英語があまり得意ではないので訳も適当ですし意味も間違っているかもしれません。その為間違いを教えている可能性も大いにあります。