- pyVRPってなに？
- pyVRPを使ってみたい…
- pyVRPで配送計画問題を解いてみたい
こんにちは！しゅんです！
この記事ではpyVRPを使って配送計画問題を解く方法を紹介したいと思います。pyVRPはオープンソースの配送計画問題（VRP）を解くソルバーです。
気になってちょっと使ってみたんですが、なんと200地点の容量制約付き配送計画問題（CVRP）を解いてみた結果、300秒で下の解が得られました。
pyVRPを日本語で紹介しているサイトが無かったので、この記事でpyVRPをpythonで動かす方法について紹介したいと思います。
それではやっていきましょう！
【Udemy講座公開のお知らせ】
このたびUdemyで数理最適化の講座を公開しました！この講座は「数理最適化を勉強してみたいけど数式が多くて難しい…」という方向けに、どうやって最適化問題を定式化すれば良いかを優しく丁寧に解説しています！
pyVRPってなに？
pyVRPは配送計画問題を解くソルバーで、pythonで利用することができます。対応している配送計画問題は
・容量制約付き配送計画問題（CVRP）
・時間枠付き配送計画問題（VRPTW）
・マルチデポ配送計画問題（MDVRP）
・同時配送集荷を考慮した配送計画問題（VRPSPD）
など色々あります。
pyVRPは厳密解を求めるソルバーではなく、短時間で高品質な解を得ることを目的としたソルバーです。またpyVRPはオープンソースでライセンスはMIT licenseです。
pythonでpyVRPを使ってCVRPを解いてみた！
それではpythonでpyVRPを使って容量制約付き配送計画問題（CVRP）を解いてみたいと思います。pyVRPのインストールはpipで行うのが一番早いと思います。
pip install pyvrpコードの全文
# --------------------------
# 必要なライブラリのインポート
# --------------------------
import random
import numpy as np
import matplotlib.pyplot as plt
from pyvrp import Model
from pyvrp.stop import MaxRuntime, MaxIterations, MultipleCriteria
from pyvrp.plotting import plot_coordinates, plot_solution
# --------------------------
# パラメータの設定
# --------------------------
num_customers = 200 # 顧客数
vehicle_count = 20 # 車両台数
vehicle_capacity = 100 # 車両の容量
max_runtime_sec = 300 # 最大実行時間（秒）
max_iterations = 100000 # 最大反復回数
np.random.seed(42) # 再現性のためのシード値設定
# --------------------------
# 座標と需要量の生成
# --------------------------
# デポ + 顧客の座標を一様分布でランダム生成（0〜1000の範囲）
coords = np.random.uniform(0, 1000, size=(num_customers + 1, 2))
# デポは中央（500, 500）に固定
coords[0] = [500, 500]
# 顧客の需要量をランダムに設定（1〜10）
demands = [0] + [random.randint(1, 10) for _ in range(num_customers)]
# --------------------------
# モデルの構築
# --------------------------
m = Model()
# 車両タイプの登録
m.add_vehicle_type(vehicle_count, capacity=vehicle_capacity)
# デポの追加
depot = m.add_depot(x=coords[0][0], y=coords[0][1])
# 顧客（クライアント）の追加
clients = [
m.add_client(x=coords[i][0], y=coords[i][1], delivery=demands[i])
for i in range(1, num_customers + 1)
]
# すべての地点（デポ + 顧客）
locations = [depot] + clients
# --------------------------
# 距離（エッジ）の定義
# --------------------------
for frm in locations:
for to in locations:
# マンハッタン距離を使用
distance = abs(frm.x - to.x) + abs(frm.y - to.y)
m.add_edge(frm, to, distance=distance)
# --------------------------
# 問題の可視化（初期状態）
# --------------------------
plt.figure(figsize=(12, 12))
ax = plt.gca()
plot_coordinates(m.data(), ax=ax)
plt.title('Vehicle Routing Problem - 200 Customers')
plt.savefig('vrp_problem.png') # ファイルとして保存
plt.close()
# --------------------------
# 停止条件の設定
# --------------------------
stop_criterion = MultipleCriteria([
MaxRuntime(max_runtime=max_runtime_sec),
MaxIterations(max_iterations=max_iterations)
])
# --------------------------
# 問題の解決
# --------------------------
print("問題を解いています...")
res = m.solve(stop=stop_criterion, display=True)
# --------------------------
# 解の可視化
# --------------------------
plt.figure(figsize=(12, 12))
ax = plt.gca()
plot_solution(res.best, m.data(), ax=ax)
plt.title('Solution - 200 Customers')
plt.savefig('vrp_solution.png') # ファイルとして保存
plt.close()
上から順番に一つずつ解説していきます。
必要なライブラリのインポート
# --------------------------
# 必要なライブラリのインポート
# --------------------------
import random
import numpy as np
import matplotlib.pyplot as plt
from pyvrp import Model
from pyvrp.stop import MaxRuntime, MaxIterations, MultipleCriteria
from pyvrp.plotting import plot_coordinates, plot_solutionこの部分では必要なライブラリのインポートを行っています。「random」、「numpy」、「matplotlib」はよく使うライブラリなので特に問題ないと思います。
下3行がpyVRPに関連するものです。これらは後々詳しく説明しますが、用途をザックリと説明すると
Model：モデルを用意するために使う
MaxRuntime、MaxIterations、MultipleCriteria：計算の停止に使う
plot_coordiantes、plot_solution：問題、結果の可視化に使う
って感じです。
パラメータの設定
# --------------------------
# パラメータの設定
# --------------------------
num_customers = 200 # 顧客数
vehicle_count = 20 # 車両台数
vehicle_capacity = 100 # 車両の容量
max_runtime_sec = 300 # 最大実行時間（秒）
max_iterations = 100000 # 最大反復回数
np.random.seed(42) # 再現性のためのシード値設定ここではパラメータの設定を行っています。ここは特に新しいことをやっているわけではないので説明は必要ないと思います。ザックリ説明すると
num_customers：顧客の数（上の例では200人に配送する）
vehicle_count：車両の数（上の例では20台の車両）
vehicle_capacity：車両の容量制限（上の例では最大100まで）
max_runtime_sec：最大計算時間（上の例では300秒で打ち切る）
また再現性を担保するために「np.random.seed(42)」でシード値を設定しています。これで毎回同じ結果が得られるようになります。
座標と需要量の生成
# --------------------------
# 座標と需要量の生成
# --------------------------
# デポ + 顧客の座標を一様分布でランダム生成（0〜1000の範囲）
coords = np.random.uniform(0, 1000, size=(num_customers + 1, 2))
# デポは中央（500, 500）に固定
coords[0] = [500, 500]
# 顧客の需要量をランダムに設定（1〜10）
demands = [0] + [random.randint(1, 10) for _ in range(num_customers)]ここでは各顧客の座標と需要量を生成しています。まず各顧客とデポの座標をランダムに生成しています。\([0,1000] \times [0,1000]\)上にランダムに生成します。ただし座標の値は整数です。
# デポ + 顧客の座標を一様分布でランダム生成（0〜1000の範囲）
coords = np.random.uniform(0, 1000, size=(num_customers + 1, 2))pyVRPは全ての数値入力を自動的に整数に変換します。そのため、データに小数値が含まれている場合は、事前にスケーリングして整数に変換する必要があります。そうしないと、予期しない動作が発生する可能性があります。
デポの座標は中央に設定します。
# デポは中央（500, 500）に固定
coords[0] = [500, 500]
また顧客の需要量は1～10の間のいずれかの整数にランダムに設定します。
# 顧客の需要量をランダムに設定（1〜10）
demands = [0] + [random.randint(1, 10) for _ in range(num_customers)]モデルの構築
# --------------------------
# モデルの構築
# --------------------------
m = Model()
# 車両タイプの登録
m.add_vehicle_type(vehicle_count, capacity=vehicle_capacity)
# デポの追加
depot = m.add_depot(x=coords[0][0], y=coords[0][1])
# 顧客（クライアント）の追加
clients = [
m.add_client(x=coords[i][0], y=coords[i][1], delivery=demands[i])
for i in range(1, num_customers + 1)
]
# すべての地点（デポ + 顧客）
locations = [depot] + clientsここではモデルの構築をしています。pyvrp専用のコードがたくさん出てきているので少し詳しく解説します。
m = Model()まずここでモデルを用意しています。イメージとしては箱を用意している感じです。ここに顧客とか車両とかの情報をどんどん入れていくみたいな感じです。
# 車両タイプの登録
m.add_vehicle_type(vehicle_count, capacity=vehicle_capacity)モデルに車両を追加するときは「.add_vehicle_type()」を使います。引数は下図のように色々ありますが、今回は一番シンプルなCVRPを解くので、「num_available」と「capacity」の2つだけを設定しています。
「num_available」は何台の車両を使って良いのか（今回は20台）を設定します。「capacity」は各車両の容量（今回は100）を設定します。
# デポの追加
depot = m.add_depot(x=coords[0][0], y=coords[0][1])デポを追加するためには「.add_depot()」を使います。引数はデポのx座標とy座標と名前です。座標は設定しなければエラーが出ますが、名前はデフォルトで空白になっているので別に設定しても設定しなくても良いです。
# 顧客（クライアント）の追加
clients = [
m.add_client(x=coords[i][0], y=coords[i][1], delivery=demands[i])
for i in range(1, num_customers + 1)
]荷物を届ける顧客を追加するためには「.add_client()」を使います。顧客は何人もいるのでリストにしました。リストの各要素が顧客を表しています。（デポは含まない）
引数は非常にたくさんありますが、今回はシンプルなCVRPを解くのでx座標とy座標と荷物の需要量の3つだけ設定します。
上の画像を見ると分かるように、荷物の集荷、サービスにかかる時間、時間枠など色々な条件を設定できます。
# すべての地点（デポ + 顧客）
locations = [depot] + clientsここではデポと顧客をまとめて「location」という1つのリストを作成しています。ここはpyVRPと関係ないので説明は省略します。
距離（エッジ）の定義
# --------------------------
# 距離（エッジ）の定義
# --------------------------
for frm in locations:
for to in locations:
# マンハッタン距離を使用
distance = abs(frm.x - to.x) + abs(frm.y - to.y)
m.add_edge(frm, to, distance=distance)
ここでは2地点間の距離（エッジ）を定義しています。for文で「location」から「frm」と「to」を取ってきて、それらのx座標とy座標からマンハッタン距離を計算しています。
マンハッタン距離にしているのは、距離を整数にするためです。前述の通りpyVRPは全ての数値入力を自動的に整数に変換します。そのため、データに小数値が含まれている場合は、事前にスケーリングして整数に変換する必要があります。そうしないと、予期しない動作が発生する可能性があります。
location中の各要素はデポもしくは顧客です。「add_depot()」の中身を見ると分かるように、デポのクラスは「Depot」です。「Depot」クラスは3つの属性（attribute）を持っていて、それぞれ「x」「y」「name」です。（これは「add_depot()」メソッドの引数と一緒）従って「Depot」クラスのインスタンス「frm」のx座標にアクセスしたい場合は「frm.x」と書けばよいという訳です。
顧客の場合も同様に「Client」クラスの属性に「x」、「y」があるので、「Client」クラスのインスタンス「frm」のx座標にアクセスしたい場合も「frm.x」と書けば良いです。
問題の可視化
# --------------------------
# 問題の可視化（初期状態）
# --------------------------
plt.figure(figsize=(12, 12))
ax = plt.gca()
plot_coordinates(m.data(), ax=ax)
plt.title('Vehicle Routing Problem - 200 Customers')
plt.savefig('vrp_problem.png') # ファイルとして保存
plt.close()ここでは問題の可視化を行っています。可視化はmatplotlibで行うんですが、pyVRPはデフォルトで問題を図示してくれるメソッドの「plot_coordinates」が備わっています。なので自分たちはこれを呼び出すだけで下図のように問題を可視化できます。
plot_coordinateメソッドの中身は普通にmatplotlibでコーディングされています。
停止条件の設定
# --------------------------
# 停止条件の設定
# --------------------------
stop_criterion = MultipleCriteria([
MaxRuntime(max_runtime=max_runtime_sec),
MaxIterations(max_iterations=max_iterations)
])ここではモデルの計算を停止する条件を設定しています。上のコードでは停止条件として「計算時間」と「反復回数」の2つを設定しています。
最大計算時間を設定したいときは「MaxRuntime()」クラスを使用します。
最大反復回数を設定したいときは「MaxIterations()」クラスを使用します。
そして今回のように停止条件が複数あるときは「MultipleCriteria()」クラスを使って1つにまとめます。
問題の解決
# --------------------------
# 問題の解決
# --------------------------
print("問題を解いています...")
res = m.solve(stop=stop_criterion, display=True)ここでは問題を解いています。問題を解くときは「Model()」クラスの「.solve()」メソッドを使います。引数の「stop」は計算の停止条件を、「display」はログをコンソールに表示するかどうかをそれぞれ表しています。
ログは下のように表示されます。
解の可視化
# --------------------------
# 解の可視化
# --------------------------
plt.figure(figsize=(12, 12))
ax = plt.gca()
plot_solution(res.best, m.data(), ax=ax)
plt.title('Solution - 200 Customers')
plt.savefig('vrp_solution.png') # ファイルとして保存
plt.close()ここでは得られた解の可視化をしています。問題の可視化の時と同じようにmatplotlibで行うんですが、pyVRPはデフォルトで解を図示してくれるメソッドの「plot_solutions」が備わっています。なので自分たちはこれを呼び出すだけで下図のように問題を可視化できます。
おわりに
今回はpyVRPを使って配送計画問題を解く方法を紹介しました！
今後もこのような組合せ最適化に関する記事を書いていきます！
最後までこの記事を読んでくれてありがとうございました。
普段は組合せ最適化の記事を書いてたりします。
ぜひ他の記事も読んでみてください！
このブログの簡単な紹介はこちらに書いてあります。
興味があったら見てみてください。
このブログでは経営工学を勉強している現役理系大学生が、経営工学に関することを色々話していきます！
ぼくが経営工学を勉強している中で感じたことや、興味深かったことを皆さんと共有出来たら良いなと思っています。
そもそも経営工学とは何なのでしょうか。Wikipediaによると
経営工学（けいえいこうがく、英: engineering management）は、人・材料・装置・情報・エネルギーを総合したシステムの設計・改善・確立に関する活動である。そのシステムから得られる結果を明示し、予測し、評価するために、工学的な分析・設計の原理・方法とともに、数学、物理および社会科学の専門知識と経験を利用する。
引用元 : 経営工学 – Wikipedia
長々と書いてありますが、要は経営、経済の課題を理系的な観点から解決する学問です。
