【作業分析3】トラックドックでの待機時間や待機台数をSimpyでシミュレーション

2023年10月16日

前回はトラックの平均到着間隔や1台当たりの平均荷卸し時間がわかっている場合に、トラックドックでの待機時間が平均何分になるかをExcelPythonを使ってシミュレーションしました。

【作業分析1】M/M/1モデルとM/G/1モデルの待ち時間をExcelでシミュレーション

【作業分析2】M/M/1モデルとM/G/1モデルの待ち時間をPythonでシミュレーション

 

ただこれでは待ち台数のシミュレーションができません。

いやできないことはないのですが、トリッキーなロジックをコーディングしないとだめで面倒です。

そこで、離散事象シミュレーションのライブラリーとして有名なSimpyを使って待ち行列のシミュレーションを行ってみました。

これなら各トラックが何分待ったかだけでなく、行列に並んだ時点で前に何台あったかまでシミュレーションすることができます。

 

Simpyとは?

物流現場では貨物は複数の処理がなされます。

例えば、入荷された貨物は検査された後、保管ロケーションに移動され、何日か後にはピッキングされます。

作業者も同じで、入荷貨物を検査したり、保管ロケーションへ移動したり、ピッキング作業を行ったりします。

つまりモノもヒトも時間的に不連続に、いろいろな処理をされたり処理したりします。

このようなのを離散事象システムといいます。

Simpyは離散事象システムのシミュレーションに使えるPythonのライブラリーです。

 

決定論的モデルと確率論的モデル

システムには決定論的モデルと確率論的モデルがあります

何回シミュレーションしても結果が同じモデルを決定論的、シミュレーションするたびに結果が変わるモデルを確率論的と言ったりします。

物流現場では同じ作業者が作業していても、毎日貨物が来るタイミングや物量も異なるので、作業が終わる時間は通常異なります。

つまり確率論的モデルです。

このように現実世界を近似したモデルは通常確率論的モデルになりますが、Simpyの使い方に少しずつ慣れていくために、まずは決定論的モデルのシミュレーションをしてみます。

 

Simpyで決定論的モデルのシミュレーション

トラックの荷卸しをする場面を考えます。

トラックバースは1つだけで、トラックバースが空になったらピッタリ15分後に貨物を積んだトラックが入ってきます。

荷卸しは15分間でピッタリ1.25台分できます。

現実的ではありませんが、このようなケースを考えます。

 

まずSimpyとNumpyをインポートします。

import simpy
import numpy  as np

次にトラックの到着間隔と荷卸しの時間間隔を定義しますが、ここでは15分間を単位時間とします。

するとトラックの到着間隔は単位時間当たり1台来るので、

1単位時間/1台=1

です。

同様に荷卸しの時間間隔は1単位時間/1.25台=0.です。

またシミュレーション時間は10,000単位時間として、これもSIM_TIMEとして定義しておきます。

# パラメータ設定
ARRIVAL_MEAN  = 1/1
SERVICE_MEAN  = 1/1.25
SIM_TIME  = 10000

次に、トラックが到着してから荷卸しが終わって出発する時間までの時間を記録しておくために、これらの結果を格納するための空のリストを作っておきましょう。

# リスト初期化
arrivals, departures = [], []
duration_in_system  = []

ここからSimpy 特有のコードを書いていきます。

Simpyには型があって、下記の3行は最低限必要なコードのようです。

# シミュレーション環境を作る
env = simpy.Environment()
# プロセスを定義する
env.process(operation(env))
# シミュレーションを実行する
env.run(until = SIM_TIME)

その上でプロセスを関数として定義します。

ここではoperationという関数名で、各工程のフローを書いていきます。

これにも型があるようで、下記のように書きます。

while True:
        時間の記録
        yield env.timeout(工程1で費やす時間)
        時間の記録
        yield env.timeout(工程2で費やす時間)
        時間の記録

これで、工程1が始まった時間工程1が終わった時間(=工程2が始まった時間)、工程2が終わった時間が記録されます。

そして先に定義した10,000単位時間になるまでシミュレーションが実行されます。

 

今回はトラックの到着を工程1、荷卸しを工程2として次のようなコードが書けます。

def operation(env):
    truck_number = 0
    while True:
        yield env.timeout(ARRIVAL_MEAN) # 前のトラックの荷卸しが終わってから次のトラックが到着するまでの時間
        time_of_arrival = env.now
        arrivals.append(time_of_arrival)
        truck_number += 1
        print("{:.0f} arrived at {:.2f}".format(truck_number, env.now))

        yield env.timeout(SERVICE_MEAN) # 荷卸しにかかる時間
        time_of_departure = env.now
        departures.append(time_of_departure)
        print("{:.0f} finished unloading at {:.2f}".format(truck_number, env.now))

以上のコードを実行すると、次のような結果が表示されます。

1 arrived at 1.00
1 finished unloading at 1.80
2 arrived at 2.80
2 finished unloading at 3.60
3 arrived at 4.60
3 finished unloading at 5.40
4 arrived at 6.40
…
5554 finished unloading at 9997.20
5555 arrived at 9998.20
5555 finished unloading at 9999.00
5556 arrived at 10000.00

1台目のトラックは1単位時間後に到着し、1.8単位時間後に荷卸しが終了しました。

2台目のトラックは2.8単位時間後に到着し、3.6単位時間後に荷卸しが終了しました。

これが延々と繰り返され、10,000単位時間後に5,556台目のトラックが到着したところでシミュレーションが終了しています。

尚、この結果は何回シミュレーションしても同じになります。

 

Simpyで確率論的モデルのシミュレーション

次に、トラックの到着間隔や荷卸しの時間間隔を確率変数にしたシミュレーションを実行してみましょう。

毎回スピードが確率的に変わる場合のシミュレーションということです。

 

到着時間と作業時間を確率的にする

そのために乱数を使いますが、全くデタラメの値の乱数ではなく、平均はわかっているものとします。

平均がわかっている場合のランダム乱数を生成する時には、指数分布を使います。

従って、工程フローを定義するoperation関数の中身は次のようになります。

def operation(env):
    truck_number = 0
    while True:
        yield env.timeout(random.expovariate(ARRIVAL_MEAN)) # トラックの到着間隔を指数分布乱数にする
        time_of_arrival = env.now
        arrivals.append(time_of_arrival)
        truck_number += 1
        print("{:.0f} arrived at {:.2f}".format(truck_number, env.now))

        yield env.timeout(random.expovariate(SERVICE_MEAN)) # 荷卸しにかかる時間を指数分布乱数にする
        time_of_departure = env.now
        departures.append(time_of_departure)
        print("{:.0f} finished unloading at {:.2f}".format(truck_number, env.now))

 結果は次のようになりました。

1 arrived at 0.78
1 finished unloading at 1.33
2 arrived at 1.51
2 finished unloading at 1.86
3 arrived at 3.51
3 finished unloading at 4.15
4 arrived at 4.24
…
5584 arrived at 9997.37
5584 finished unloading at 9997.49
5585 arrived at 9998.54
5585 finished unloading at 9998.73

 今度はトラックの到着間隔や荷卸し時間が毎回変わっています。

また、シミュレーションをするたびに、この結果は変わります

 

リソースを定義して待ち時間を考慮する

次に待ち時間をシミュレーションできるようにしてみましょう。

待ち時間が発生するのは、リソースが限られているためです。

Simpyではリソースクラスを使って、この限られたリソースを定義できます。

今回の例ではトラックドックがリソースになり、キャパシティはドックが1つしかないので1です。

まずはキャパシティを1として、リソースクラスをインスタンス化しておきます。

変数ドックの数を表すNUMBER_OF_BAYSには、予め1を代入しておきます。

# リソースを定義する
operation_capa = simpy.Resource(env, capacity = NUMBER_OF_BAYS)

そしてこのインスタンスを使って、プロセスを記述する関数の中で次のように記述します。

with インスタンス名.request() as req:
    行列に並んだ時間の記録
    length = len(インスタンス名.queue) # 行列に並んだ時の行列の長さを取得
    yield req
    行列から出た時間の記録

こうすることにより、行列に並び始めた時間、その時の行列の長さ、行列から出た時の時間を取得することができます。

今回の例では次のようにコーディングしました。

def operation(env, operation_capa):
    truck_number = 0
    while True:
        yield env.timeout(random.expovariate(ARRIVAL_MEAN))
        time_of_arrival = env.now
        arrivals.append(time_of_arrival)
        truck_number += 1
        print("{:.0f} arrived at {:.2f}".format(truck_number, env.now))

        with operation_capa.request() as req:
            print("{:.0f} entered queue at {:.2f}".format(truck_number, env.now))
            length = len(operation_capa.queue)
            len_in_queue.append(length)
            queue_in = env.now
            tm_in_queue.append(queue_in)

            yield req
            print("{:.0f} left queue at {:.2f}".format(truck_number, env.now))
            queue_out = env.now
            tm_out_queue.append(queue_out)

            yield env.timeout(random.expovariate(SERVICE_MEAN))
            time_of_departure = env.now
            departures.append(time_of_departure)
            print("{:.0f} finished unloading at {:.2f}".format(truck_number, env.now))

結果は次のようになりました。

1 arrived at 1.02
1 entered queue at 1.02
1 left queue at 1.02
1 finished unloading at 1.11
2 arrived at 2.65
2 entered queue at 2.65
2 left queue at 2.65
2 finished unloading at 3.40
3 arrived at 3.41
…
5570 finished unloading at 9994.87
5571 arrived at 9995.78
5571 entered queue at 9995.78
5571 left queue at 9995.78
5571 finished unloading at 9997.57

各トラックについての到着時間、行列INした時間、行列OUTした時間、荷卸しが終わった時間が記録されていますが、何かが違います。

そうです、待ち時間がすべてゼロなのです。

これは前のトラックの荷卸しが終了しないと次のトラックが到着しないからです。

実際には、前のトラックの荷卸しに時間がかかっていれば、次のトラック、いや時には次の次のトラックまでもが容赦なく到着してくるはずです。

 

M/M/1モデルをシミュレーションする

それをシミュレーションするためには、到着のプロセスを独立させて10,000単位時間になるまで延々と繰り返します。

そして、トラックドックに並んで荷卸しが終わるまでのプロセスを別に作り、トラックが到着したらそちらへ受け渡すようにしてやります。

そうすればトラックが容赦なく到着する環境で、ドックが空いていればすぐに荷卸し、空いていなければ空くまで待ってから荷卸しするという状況を再現できます。

具体的にはメインのプログラムでtruck_arrivalというプロセス(トラックが到着するプロセス)を定義して10,000単位時間になるまで実行し、truck_arrivalの中で別のoperationというプロセス(ドックが空くまで待ってから荷卸しをするプロセス)を定義して受け渡してやることにします。

# パラメータ設定
NUMBER_OF_BAYS = 1
ARRIVAL_MEAN  = 1
SERVICE_MEAN  = 1.25
SIM_TIME  = 10000

# リスト初期化
arrivals, departures = [], []
duration_in_queue, duration_in_service, duration_in_system  = [], [], []
tm_in_queue, tm_out_queue = [], []
len_in_queue = []

def truck_arrival(env, operation_capa):
    truck_number = 0
    while True:
        yield env.timeout(random.expovariate(ARRIVAL_MEAN))
        time_of_arrival = env.now
        arrivals.append(time_of_arrival)
        truck_number += 1
        print("{:.0f} arrived at {:.2f}".format(truck_number, env.now))
        env.process(operation(env, operation_capa, truck_number, time_of_arrival))

def operation(env, operation_capa, truck_number, time_of_arrival):
    with operation_capa.request() as req:
        print("{:.0f} entered queue at {:.2f}".format(truck_number, env.now))
        length = len(operation_capa.queue)
        len_in_queue.append(length)
        queue_in = env.now
        tm_in_queue.append(queue_in)

        yield req
        print("{:.0f} left queue at {:.2f}".format(truck_number, env.now))
        queue_out = env.now
        tm_out_queue.append(queue_out)

        yield env.timeout(random.expovariate(SERVICE_MEAN))
        time_of_departure = env.now
        departures.append(time_of_departure)
        print("{:.0f} finished unloading at {:.2f}".format(truck_number, env.now))

        time_in_system = time_of_departure - time_of_arrival
        duration_in_system.append(time_in_system)
        time_in_service = time_of_departure - queue_out
        duration_in_service.append(time_in_service)
        time_in_queue = queue_out - queue_in
        duration_in_queue.append(time_in_queue)

# シミュレーション環境を作る
env = simpy.Environment()
# リソースを定義する
operation_capa = simpy.Resource(env, capacity = NUMBER_OF_BAYS)
# プロセスを定義する
env.process(truck_arrival(env, operation_capa))
# シミュレーションを実行する
env.run(until = SIM_TIME)

結果は次のようになりました。

1 arrived at 0.81
1 entered queue at 0.81
1 left queue at 0.81
2 arrived at 1.41
2 entered queue at 1.41
1 finished unloading at 1.56
2 left queue at 1.56
2 finished unloading at 1.87
3 arrived at 6.99
3 entered queue at 6.99
3 left queue at 6.99
4 arrived at 8.34
4 entered queue at 8.34
3 finished unloading at 8.43
…

1台目が荷卸しをしている最中の1.41単位時間に2台目のトラックが到着し、ドックが空く1.56単位時間まで待っている様子がわかります。

3台目は2台目の荷卸しが終了してから6.99単位時間に到着していますが、4台目はまたもや3台目の荷卸し中に到着しています。

このようにトラックがドックで待っている様子がシミュレーションできていることがわかります。

 

先のコードではduration_in_queueのリストに行列での待ち時間、duration_in_serviceに荷卸し時間、duration_in_systemにその合計(各トラック当たりの所要時間)を格納していますので、それらの合計を計算してみましょう。

np.mean(duration_in_queue), np.mean(duration_in_service), np.mean(duration_in_system)
>>
(3.3755474385469806, 0.8009464850021448, 4.176493923549126)

これら3つの値の理論値はそれぞれ3.20.84.0なので、それらに近しい値になっています。

【待ち行列|M/M/1モデル】5つの定理の証明方法と物流への適用事例

また各トラックがドックに入るまでに前に何台のトラックが並んでいたかは、len_in_queueのリストに格納されています。

その平均と最大値を求めてみたら、次のようになりました。

np.mean(len_in_queue), np.max(len_in_queue)
>>
(4.192713989534998, 29)

平均4台が前に並んでいたけれども、多い時には29台も並んでいたことがわかりました。

M/M/1モデルの公式では平均4台というのは分かれど、最大29台というのは分かりません。

シミュレーションの有用性がここにあります。

次回はこの結果を動画で表示させてみましょう。

【作業分析4】トラックの待機台数のグラフをPythonでアニメーション(1)