📚

PuLPとOR-Tools (CP-SAT) の実装を見比べながらOR-Toolsの入門

に公開

この記事は数理最適化 Advent Calendar 2025の 16 日目の記事です。

はじめに

結構前に仕事の方でortoolsを使ったことがあったんですが、公式ドキュメント以外、案外実装例がなかったこともあって、まあそこそこ困りました。
pulpはめっちゃ出てくるんですけどね。
今はChatGPTに聞けば出てくるかもですが、備忘録も兼ねて両者の書き方の違いをまとめておきます。

ちなみにドキュメントはこれ: https://or-tools.github.io/docs/pdoc/ortools/sat/python/cp_model.html
これも地味に探しづらいなあーとか勝手に思ってます。

1. 基本的な違い

まず、基本的な違いから見ていきます。PuLPは混合整数計画法(MIP)、OR-ToolsのCP-SATは制約プログラミング(CP)というアプローチです。

PuLP(混合整数計画法)

PuLPは数学的な式をほぼそのまま書けるのが特徴です。線形計画法や混合整数計画法の教科書に出てくる式をそのままコードに落とせる感じですね。

import pulp

# 問題定義
problem = pulp.LpProblem("sample", pulp.LpMinimize)

# 変数定義
x = pulp.LpVariable("x", lowBound=0, upBound=10, cat='Integer')
y = pulp.LpVariable("y", cat='Binary')

# 目的関数
problem += 2*x + 3*y

# 制約
problem += x + y <= 5
problem += 2*x - y >= 1

# 求解
problem.solve()

変数を定義したら、+=で制約を追加していくだけです。数式を書く感覚でモデリングできます。

OR-Tools CP-SAT(制約プログラミング)

一方、OR-Toolsは少し違ったアプローチです。変数の型ごとにメソッドが用意されていて、制約も専用のメソッドで追加します。

from ortools.sat.python import cp_model as cp

# モデル定義
model = cp.CpModel()

# 変数定義
x = model.NewIntVar(0, 10, "x")
y = model.NewBoolVar("y")

# 制約を追加
model.Add(x + y <= 5)
model.Add(2*x - y >= 1)

# 目的関数
model.Minimize(2*x + 3*y)

# ソルバー実行
solver = cp.CpSolver()
status = solver.Solve(model)

# 結果取得
if status == cp.OPTIMAL:
    print(f"x = {solver.Value(x)}")
    print(f"y = {solver.Value(y)}")

見た目はそこまで変わらないですが、変数定義がNewIntVarNewBoolVarのようにメソッド呼び出しになっています。制約もAdd()メソッドを使います。

慣れるまではpulpの方が直感的で書きやすかなと。
pulpは目的関数のために説明変数を操作していく感じ。
cpは説明変数がどうのってより、制約を一つずつ追加していく感じ。
とはいえ、そんなに大きな差はないです。

2. 変数定義の比較

次に、変数定義の書き方を見ていきます。基本的にはPuLPもOR-Toolsも似たような感じですが、細かい記法が違います。

整数変数

整数変数の定義はこんな感じです。

PuLP:

x = pulp.LpVariable("x", lowBound=0, upBound=100, cat='Integer')

OR-Tools:

x = model.NewIntVar(0, 100, "x")

PuLPはcat='Integer'で指定するのに対し、OR-ToolsはNewIntVarというメソッド名自体が型を表しています。引数の順番も違って、OR-Toolsは(下限、上限、変数名)の順です。

バイナリ変数

バイナリ変数(0か1の変数)も同様です。

PuLP:

y = pulp.LpVariable("y", cat='Binary')

OR-Tools:

y = model.NewBoolVar("y")

OR-ToolsではNewBoolVarという専用メソッドがあります。バイナリ変数をよく使うので、専用メソッドがあるのは便利ですね。

連続変数

PuLP:

z = pulp.LpVariable("z", lowBound=0, upBound=1, cat='Continuous')

OR-Tools:

# CP-SATは整数変数のみ扱う
# 連続変数は適切にスケーリングして整数化する必要がある
# 例:0.01刻みで0〜1の範囲なら
z = model.NewIntVar(0, 100, "z_scaled")  # 実際の値は z/100

OR-ToolsのCP-SATは整数変数専用なので、連続変数が必要な場合はスケーリングが必要になります。例えば、0.01刻みで扱いたいなら100倍して整数変数にして、後で100で割り戻すという感じです。この辺はちょっと面倒ですね。

自分が実装した車両スケジューリングでは、時間を分単位で扱っていたので、この問題には当たりませんでした。もし分単位や小数点を含む時間を扱う場合は、秒単位に変換してから扱うといいかと思います。

3. 制約の書き方

基本的な線形制約

次に制約の書き方です。基本的な線形制約はPuLPもOR-Toolsもほぼ同じように書けます。

PuLP:

# 数式をそのまま書ける
problem += 2*x + 3*y <= 10
problem += x - y == 0
problem += x >= 5

OR-Tools:

# Add()メソッドを使う
model.Add(2*x + 3*y <= 10)
model.Add(x - y == 0)
model.Add(x >= 5)

この辺の基本的な線形制約は、PuLPとほぼ同じ感覚で書けますね。+=Add()かの違いだけで、不等式の書き方自体は同じです。

sumの扱い

複数の変数の合計を扱う場合も、基本的には同じです。

PuLP:

variables = [pulp.LpVariable(f"x_{i}", cat='Binary') for i in range(10)]

problem += sum(variables) == 5

OR-Tools:

variables = [model.NewBoolVar(f"x_{i}") for i in range(10)]

model.Add(sum(variables) == 5)
# or
model.Add(cp.LinearExpr.Sum(variables) == 5)

OR-ToolsでもPythonのsum()が使えます。ただ、大規模な問題の場合はLinearExpr.Sumを使った方が効率的らしいです。
自分はこれが原因で、なんか定式化にすげー時間かかってんなっていう場面に出会しました。
あとはまあ、せっかくの型定義がsumだとされなくなるので、ライブラリ付属を使った方が何かとお得かなとは思います。
とはいえ、面倒とか普通の規模ならsum()で十分だと思います。

ExactlyOne制約(排他的選択)

「複数の選択肢から必ず1つだけ選ぶ」という制約は、スケジューリング問題でよく出てきます。

PuLP:

problem += sum(variables) == 1

OR-Tools:


model.AddExactlyOne(variables)

# or
model.AddAtMostOne(variables)   # 最大1つ

# or
model.AddAtLeastOne(variables)  # 最低1つ

OR-Toolsには論理制約のための専用メソッドが用意されていて、単にsum()==1と書くよりソルバーが効率的に処理してくれます。

だし、より素直な表現になっているので、昨今のLLM時代には相性がいいのかもしれませんね。

4. OR-Toolsの特有機能

orotoolsでは制約を簡単に書ける機能があります。

4.1 条件付き制約(OnlyEnforceIf)

PuLPで条件付き制約を書こうとすると、ビッグMというテクニックを使う必要があります。
ちゃんとやればMの見積もりはキチンでとできることがおおいし、ちゃんとやらないとパフォーマンスに影響が出たりします。

PuLP:

x = pulp.LpVariable("x", lowBound=0, upBound=100, cat='Integer')
y = pulp.LpVariable("y", cat='Binary')

M = 1000  # 十分大きな数

# y==1のとき x >= 50、y==0のとき制約なし
problem += x >= 50 - M * (1 - y)
problem += x <= 50 + M * (1 - y)

とはいえ、ちょこっとやるのにこんなの考えてられるかって人も多いかもです。
OR-Toolsではもっと直観的ににかけたりします。

OR-Tools(OnlyEnforceIf):

x = model.NewIntVar(0, 100, "x")
y = model.NewBoolVar("y")

# y==1(True)のときのみ x >= 50
model.Add(x >= 50).OnlyEnforceIf(y)

# y==0(False)のときのみ x < 50
model.Add(x < 50).OnlyEnforceIf(y.Not())

Mが消えてますね。
次の節でも出てきますが、Mを消せるような場面が多く出てきます。
特に大きな問題を解く時はLP緩和なんか考えたりすると思うんですが、こう言ったグローバルな定数が現れるとすんごいめんどいです。
この辺りはヒューリュスティックがいきる場面かなと思います。

4.2 Interval変数とNoOverlap制約

このあたりが便利メソッドというか、一番の違いかなと。
スケジューリング問題でよくある「複数のタスクを1つのリソースで実行するけど、タスクは重複できない」という制約を見てみます。

PuLPで実装するとこうなります。

PuLP:

# タスク i の開始時間
start_times = [
	pulp.LpVariable(f"start_{i}", lowBound=0, upBound=100, cat='Integer')
    for i in range(n_tasks)
]

# タスクの処理時間
durations = [10, 20, 15, 25]

# タスク i が タスク j より前に実行されるかのバイナリ変数
M = 1000
for i in range(n_tasks):
    for j in range(i+1, n_tasks):
        y = pulp.LpVariable(f"order_{i}_{j}", cat='Binary')

        # y==1: タスクiがタスクjより前
        problem += start_times[i] + durations[i] <= start_times[j] + M * (1 - y)
        # y==0: タスクjがタスクiより前
        problem += start_times[j] + durations[j] <= start_times[i] + M * y

タスクが4つなら、6個のバイナリ変数と12個の制約が必要になります。タスクが10個なら45個の変数と90個の制約...という感じで、タスク数が増えると爆発的に増えます。
ここでもMが出てきます。

OR-Toolsだと以下になります。

OR-Tools(Interval変数とNoOverlap):

# タスクの開始時間
start_times = [model.NewIntVar(0, 100, f"start_{i}") for i in range(n_tasks)]

# タスクの処理時間
durations = [10, 20, 15, 25]

# Interval変数を作成
intervals = [
    model.NewIntervalVar(
	    start_times[i],
	    durations[i], 
	    start_times[i] + durations[i], 
	    f"interval_{i}"
	)
    for i in range(n_tasks)
]

# NoOverlap制約を追加
model.AddNoOverlap(intervals)

やってることはInterval変数を作って、AddNoOverlapを呼ぶだけ。

O(n²)の変数と制約が必要になったりしますが、OR-ToolsならInterval変数とNoOverlap制約だけで済みます。コードの可読性も上がってます。

例えば車両スケジューリングでは、ターミナルのリソース競合を表現するのにこのNoOverlap制約をつかうなどですね。
複数の車両が同じターミナルを使うけど、同時には使えないのような。

4.3 Optional Interval変数

さらに、「条件によってタスクが実行されるかどうか変わる」場合もあります。例えば「このタスクは必要に応じて実行する」みたいなケースです。

OR-Tools:

# タスク実行の有無を決めるバイナリ変数
is_active = model.NewBoolVar("is_active")

# Optional Interval: is_activeがTrueの時のみ有効
interval = model.NewOptionalIntervalVar(
    start=start_time,
    size=duration,
    end=start_time + duration,
    is_present=is_active,
    name="optional_interval"
)

# NoOverlapにOptional Intervalも含められる
model.AddNoOverlap([interval, other_interval])

個人的にはこれが一番便利やんとなったやつです。
is_activeという2値変数の如何によって、interval変数を使うかどうかを決められます。

実用的には、最適化を考える上で、車両がこの道路を通るとか通らないとかでx_{e} \in \{0,1\}の変数があったとして、この値をそのままis_presentへ利用できます。
通る場合は、この制約追加する。通らない場合はこの制約を追加するのような操作です。

あとは、オートマトンの制約とかまた結構使いこなしたいなあというのは結構あります。

5. ソルバーの実行と結果取得

最後に、ソルバーの実行方法と結果の取得方法を見ておきますがここはほぼ同じです。
中身だけ書いておきます。

PuLP

status = problem.solve(pulp.PULP_CBC_CMD(msg=1))

# 結果確認
if status == pulp.LpStatusOptimal:
    print(f"最適値: {pulp.value(problem.objective)}")
    for v in problem.variables():
        print(f"{v.name} = {v.varValue}")
else:
    print(f"ステータス: {pulp.LpStatus[status]}")

solve()を呼ぶだけです。デフォルトのソルバー(CBC)が使われますが、Gurobiなど他のソルバーを指定することもできます。

OR-Tools

OR-Toolsは少し詳しい情報が取れます。

solver = cp.CpSolver()

# パラメータ設定(必要に応じて)
solver.parameters.max_time_in_seconds = 300.0
solver.parameters.num_search_workers = 8  # 並列数

# 実行
status = solver.Solve(model)

# 結果確認
if status == cp.OPTIMAL:
    print(f"最適値: {solver.ObjectiveValue()}")
    print(f"x = {solver.Value(x)}")
    print(f"y = {solver.Value(y)}")

    # 統計情報
    print(f"解法時間: {solver.WallTime()}秒")
    print(f"探索ノード数: {solver.NumBranches()}")

elif status == cp.FEASIBLE:
    print("実行可能解が見つかりました(最適性は未保証)")
    print(f"現在の目的関数値: {solver.ObjectiveValue()}")

else:
    print(f"解が見つかりませんでした: {solver.StatusName()}")

6. 全体実装例:容量制約VRP

容量制約VRPを例に、PuLPとOR-Toolsの実装を比較してみます。

問題設定

配送会社が複数の顧客に荷物を配送する問題です。

前提条件:

  • デポ(倉庫)から複数の車両が出発
  • 各車両には積載容量の制限がある
  • すべての顧客を訪問し、デポに戻る

目的: 総移動距離(または時間)を最小化

主な制約:

  1. 容量制約: 車両の積載量を超えない
  2. 訪問制約: すべての顧客を1回ずつ訪問
  3. 経路制約: デポから出発してデポに戻る

import math

# 問題データ
n_customers = 10  # 顧客数
n_vehicles = 3    # 車両数
vehicle_capacity = 100  # 車両容量

# 顧客データ
customers = {
    0: {'demand': 0, 'x': 0, 'y': 0},      # デポ
    1: {'demand': 20, 'x': 10, 'y': 5},
    2: {'demand': 30, 'x': 5, 'y': 10},
    3: {'demand': 15, 'x': 15, 'y': 8},
    4: {'demand': 25, 'x': 8, 'y': 15},
    5: {'demand': 10, 'x': 20, 'y': 12},
    6: {'demand': 35, 'x': 12, 'y': 20},
    7: {'demand': 18, 'x': 3, 'y': 18},
    8: {'demand': 22, 'x': 18, 'y': 3},
    9: {'demand': 28, 'x': 7, 'y': 22},
}

# 距離計算
def calc_distance(i, j):
    return int(math.sqrt((customers[i]['x'] - customers[j]['x'])**2 +
                         (customers[i]['y'] - customers[j]['y'])**2))

distances = {(i, j): calc_distance(i, j)
             for i in customers for j in customers if i != j}

PuLP実装

まず、PuLPでの実装です。基本的なVRPの定式化です。

pulp実装
import pulp

# 問題定義
problem = pulp.LpProblem("CVRP", pulp.LpMinimize)

# 変数定義
# x[i,j,k] = 車両kが顧客iから顧客jへ移動するか
x = {}
for k in range(n_vehicles):
    for i in customers:
        for j in customers:
            if i != j:
                x[i,j,k] = pulp.LpVariable(f"x_{i}_{j}_{k}", cat='Binary')

# u[i,k] = 車両kが顧客iを出発する時点での積載量(MTZ制約用)
u = {}
for k in range(n_vehicles):
    for i in customers:
        u[i,k] = pulp.LpVariable(f"u_{i}_{k}",
                                 lowBound=0, upBound=vehicle_capacity,
                                 cat='Integer')

# 目的関数:総移動距離の最小化
problem += pulp.lpSum(distances[i,j] * x[i,j,k]
                      for k in range(n_vehicles)
                      for i in customers for j in customers if i != j)

# 制約1:各顧客は必ず1回訪問される
for i in customers:
    if i == 0:  # デポは除く
        continue
    problem += pulp.lpSum(x[i,j,k]
                         for k in range(n_vehicles)
                         for j in customers if i != j) == 1

# 制約2:各車両はデポから出発
for k in range(n_vehicles):
    problem += pulp.lpSum(x[0,j,k] for j in customers if j != 0) == 1

# 制約3:各車両はデポに戻る
for k in range(n_vehicles):
    problem += pulp.lpSum(x[i,0,k] for i in customers if i != 0) == 1

# 制約4:流量保存
for k in range(n_vehicles):
    for i in customers:
        if i == 0:
            continue
        problem += (pulp.lpSum(x[i,j,k] for j in customers if i != j) ==
                   pulp.lpSum(x[j,i,k] for j in customers if i != j))

# 制約5:容量制約とMTZ制約
M_capacity = vehicle_capacity
for k in range(n_vehicles):
    # デポでの積載量は0
    problem += u[0,k] == 0

    # デポから顧客への移動:初期積載量を設定
    for j in customers:
        if j == 0:
            continue
        # x[0,j,k] == 1 のとき、u[j,k] >= demand[j]
        problem += u[j,k] >= customers[j]['demand'] - M_capacity * (1 - x[0,j,k])

    # 顧客間の移動:積載量を累積
    for i in customers:
        if i == 0:
            continue
        for j in customers:
            if j == 0 or i == j:
                continue
            # x[i,j,k] == 1 のとき、u[j,k] >= u[i,k] + demand[j]
            problem += (u[j,k] >= u[i,k] + customers[j]['demand'] -
                       M_capacity * (1 - x[i,j,k]))

    # 各顧客での積載量の上限(容量制約)
    for i in customers:
        if i == 0:
            continue
        problem += u[i,k] <= vehicle_capacity

# 求解
solver = pulp.PULP_CBC_CMD(msg=1, timeLimit=300)
status = problem.solve(solver)

# 結果取得
print(f"目的関数値: {pulp.value(problem.objective)}")

多分こんな感じの実装になるかと。

OR-Toolsの実装

次に、OR-ToolsのCP-SATでの実装です。

OR-Toolsの実装
from ortools.sat.python import cp_model as cp

# モデル定義
model = cp.CpModel()

# 変数定義
# x[i,j,k] = 車両kが顧客iから顧客jへ移動するか
x = {}
for k in range(n_vehicles):
    for i in customers:
        for j in customers:
            if i != j:
                x[i,j,k] = model.NewBoolVar(f"x_{i}_{j}_{k}")

# u[i,k] = 車両kが顧客iを出発する時点での積載量(MTZ制約用)
u = {}
for k in range(n_vehicles):
    for i in customers:
        u[i,k] = model.NewIntVar(0, vehicle_capacity, f"u_{i}_{k}")

# 目的関数:総移動距離の最小化
total_distance = sum(distances[i,j] * x[i,j,k]
                    for k in range(n_vehicles)
                    for i in customers for j in customers if i != j)
model.Minimize(total_distance)

# 制約1:各顧客は必ず1回だけ訪問される
for i in customers:
    if i == 0:
        continue
    # すべての車両からのすべての入り辺の合計が1
    model.Add(sum(x[j,i,k] for k in range(n_vehicles)
                            for j in customers if i != j) == 1)

# 制約2:各車両はデポから出発してデポに戻る
for k in range(n_vehicles):
    model.Add(sum(x[0,j,k] for j in customers if j != 0) == 1)
    model.Add(sum(x[i,0,k] for i in customers if i != 0) == 1)

# 制約3:流量保存
for k in range(n_vehicles):
    for i in customers:
        if i == 0:
            continue
        outgoing = sum(x[i,j,k] for j in customers if i != j)
        incoming = sum(x[j,i,k] for j in customers if i != j)
        model.Add(outgoing == incoming)

# 制約4:容量制約とMTZ制約
for k in range(n_vehicles):
    # デポでの積載量は0
    model.Add(u[0,k] == 0)

    # デポから顧客への移動:初期積載量を設定
    for j in customers:
        if j == 0:
            continue
        # x[0,j,k] == 1 のとき、u[j,k] >= demand[j]
        model.Add(u[j,k] >= customers[j]['demand']).OnlyEnforceIf(x[0,j,k])

    # 顧客間の移動:積載量を累積
    for i in customers:
        if i == 0:
            continue
        for j in customers:
            if j == 0 or i == j:
                continue
            # x[i,j,k] == 1 のとき、u[j,k] >= u[i,k] + demand[j]
            model.Add(u[j,k] >= u[i,k] + customers[j]['demand']).OnlyEnforceIf(x[i,j,k])

    # 各顧客での積載量の上限(容量制約)
    for i in customers:
        if i == 0:
            continue
        model.Add(u[i,k] <= vehicle_capacity)


# ソルバー実行
solver = cp.CpSolver()
solver.parameters.max_time_in_seconds = 60.0
status = solver.Solve(model)

# 結果取得
print(f"目的関数値: {solver.ObjectiveValue()}")

解いた結果

ルート表示

目的関数値: 132.0

=== 制約違反チェック ===

  1. 顧客訪問回数:
    顧客1: 1回 ✓
    顧客2: 1回 ✓
    顧客3: 1回 ✓
    顧客4: 1回 ✓
    顧客5: 1回 ✓
    顧客6: 1回 ✓
    顧客7: 1回 ✓
    顧客8: 1回 ✓
    顧客9: 1回 ✓

  2. 車両容量制約:
    車両0:
    デポ: u[0,0] = 0
    顧客1: u[1,0] = 92 (需要=20) ✓
    顧客3: u[3,0] = 50 (需要=15) ✓
    顧客4: u[4,0] = 25 (需要=25) ✓
    顧客5: u[5,0] = 35 (需要=10) ✓
    顧客8: u[8,0] = 72 (需要=22) ✓
    → 最大積載量: 92 / 100 ✓
    車両1:
    デポ: u[0,1] = 0
    顧客2: u[2,1] = 30 (需要=30) ✓
    → 最大積載量: 30 / 100 ✓
    車両2:
    デポ: u[0,2] = 0
    顧客6: u[6,2] = 100 (需要=35) ✓
    顧客7: u[7,2] = 18 (需要=18) ✓
    顧客9: u[9,2] = 46 (需要=28) ✓
    → 最大積載量: 100 / 100 ✓

=== ルート ===
車両0: [0, 4, 5, 3, 8, 1, 0] (需要合計: 92)
車両1: [0, 2, 0] (需要合計: 30)
車両2: [0, 7, 9, 6, 0] (需要合計: 81)

問題が問題なので、こう書くとどっちも記述量変わらんように見えますね。。。
実際に自分で手を動かすといい感じ使いやすさは出てくるかなと。

7. まとめ

導入するにはやっぱり実装難易度ってのは重要かなと思います。
とはいえ、当たり前ですが、問題の種類によって向き不向きがあります。
厳密解が必要な場合は、モデリングツールとしてはここのpulpとか、他にamplとかgurobipyとかその辺りになるかなと。
ヒューリスティックでよければ、割と直観的な表現ができるので良し悪しあれ、そこまで数式を書かなくても実装できるのは障壁が低いのかなと。

Fusic 技術ブログ

Discussion