はじめに
筆者は前職にて2次元パッキング問題を扱う機会があり、振り返りも兼ねてこの記事を書きました。
この記事では、2次元矩形パッキング問題に対する基本的なアルゴリズムであるBottom-Left法（BL法）の実装と、No-Fit Polygon（NFP）を用いた高速化について解説します。
この記事で扱うこと
- BL法の理論と実装（O(n⁴)）
- NFP・IFRの理論と実装
- 走査線アルゴリズム（Find2D-BL）による高速化（O(n²logn)）
- 単純版とNFP版のベンチマーク比較
この記事で扱わないこと
- 多角形パッキング（別記事で扱う予定）
- メタヒューリスティクスによる解品質の改善（別記事で扱う予定）
実装コードはGitHubで公開しており、Streamlitによるインタラクティブなデモも用意しています。
- GitHub：https://github.com/Haru8-8/2d-packing-algorithms
- デモ：https://2d-packing-algorithms-r3lvpgajmvk2earrnhl3dk.streamlit.app/
2次元パッキング問題とは
2次元パッキング問題とは、与えられた複数のアイテムを重なりなく指定された領域（ビン）に配置する問題の総称です。
本記事では、その中でも2次元ストリップパッキング問題を扱います。幅が固定された矩形のビンに、複数の矩形アイテムを重なりなく配置し、使用するビンの高さを最小化することが目的です。
この問題は工場のレイアウト設計、鉄鋼・ガラス・繊維などの素材の裁断、集積回路の設計など、様々な分野に応用があります。
一見シンプルに見えますが、ほぼすべてのパッキング問題はNP困難であることが知られており、多項式時間で最適解を求めるアルゴリズムはおそらく存在しません。そのため実用上は、近似解を高速に求めるヒューリスティクスが広く用いられています。
BL法（Bottom-Left法）
アルゴリズムの概要
BL法は、空のビンから始めて矩形を1つずつ BL点（最も下で最も左の位置） に配置していく構築型ヒューリスティクスです。アルゴリズムはシンプルで、以下の手順を繰り返します。
- 配置する矩形をあらかじめ決めた順序（面積降順など）に並べる
- 各矩形について、配置可能な位置の候補を列挙する
- 候補の中でy座標が最小、同じならx座標が最小の点（BL点）に配置する
BL実行可能点とBL点
矩形をある座標に置いたとき、以下の条件を同時に満たす座標をBL実行可能点と呼びます。
- 既配置の矩形と重ならない
- ビンからはみ出さない
- その座標から下にも左にも動かせない
BL実行可能点の中でy座標が最小、同じならx座標が最小の点がBL点です。
候補点の列挙
BL実行可能点は必ず以下のいずれかの座標で決まります。
| 種別 | 下方向に接するもの | 左方向に接するもの |
|---|---|---|
| A | 容器の底辺 | 容器の左辺 |
| B | 矩形kの上辺 | 矩形jの右辺 |
| C | 矩形jの上辺 | 容器の左辺 |
| D | 容器の底辺 | 矩形jの右辺 |
これらの組み合わせはO(n²)個あり、各候補に対してO(n)の重なり判定を行うため、1矩形あたりO(n³)、全体でO(n⁴)の計算量になります。
Pythonによる実装
def bl_candidates(i, x, y, w, h):
"""BL実行可能点の候補を列挙する。"""
cand = [(0.0, 0.0)] # 容器の左下
for j in range(i):
for k in range(j):
# 矩形jの右端 × 矩形kの上端
cand.append((x[j] + w[j], y[k] + h[k]))
# 矩形kの右端 × 矩形jの上端
cand.append((x[k] + w[k], y[j] + h[j]))
# 容器左端 × 矩形jの上端
cand.append((0.0, y[j] + h[j]))
# 矩形jの右端 × 容器底辺
cand.append((x[j] + w[j], 0.0))
return cand
def is_feasible(i, p, x, y, w, h, bin_w):
"""座標pに矩形iを配置できるか判定する。"""
px, py = p
if px < 0 or px + w[i] > bin_w or py < 0:
return False
for j in range(i):
x_overlap = max(px, x[j]) < min(px + w[i], x[j] + w[j])
y_overlap = max(py, y[j]) < min(py + h[i], y[j] + h[j])
if x_overlap and y_overlap:
return False
return True
def bl_method(rects, bin_w, sort_key='area'):
"""BL法で矩形をストリップに詰める。"""
# 面積降順に並べ替え
order = sorted(range(len(rects)),
key=lambda k: rects[k][0] * rects[k][1], reverse=True)
sorted_rects = [rects[k] for k in order]
w = [r[0] for r in sorted_rects]
h = [r[1] for r in sorted_rects]
x, y = [], []
for i in range(len(sorted_rects)):
cands = bl_candidates(i, x, y, w, h)
feasible = [p for p in cands if is_feasible(i, p, x, y, w, h, bin_w)]
bl_point = min(feasible, key=lambda p: (p[1], p[0]))
x.append(bl_point[0])
y.append(bl_point[1])
return list(zip(x, y)), sorted_rects
解の例
面積降順でBL法を実行した例です。ソート順によって充填率が変わることが確認できます。
※ GitHubリポジトリのノートブック01_bottom_left.ipynbで再現できます。
より少ない計算量で同等の解を得るために
BL法の単純実装はO(n⁴)の計算量を持つため、矩形数が増えると実行時間が急激に増加します。実際に計測すると、n=100の場合で約4秒かかります。
| n | 実行時間 |
|---|---|
| 20 | 0.0100s |
| 50 | 0.2928s |
| 100 | 4.2966s |
実用的な規模の問題（n=100以上）では、このままでは厳しい状況です。
ボトルネックはBL点の探索にあります。候補点の列挙がO(n²)、各候補への重なり判定がO(n)であるため、1矩形あたりO(n³)かかっています。
ここでNo-Fit Polygon（NFP）と走査線アルゴリズムを組み合わせることで、BL点探索を1矩形あたりO(n log n)に改善できます。結果として全体の計算量はO(n²log n)となり、単純実装と全く同じ配置結果を得ながら大幅に高速化できます。
NFP（No-Fit Polygon）とは
NFPの定義
No-Fit Polygon（NFP）とは、平面上で2つの図形の重なりを判定するためのデータ構造です。
固定された図形Pに対して、移動図形Qの 参照点（左下の頂点） が取り得る座標のうち、PとQが重なりを持つような座標の集合をNFP(P, Q)と呼びます。
NFPには以下の重要な性質があります。
- Qの参照点がNFPの内部にある → PとQは重なる
- Qの参照点がNFPの境界上にある → PとQは接触（重ならない）
- Qの参照点がNFPの外部にある → PとQは重ならない
つまり、「Qの参照点がNFPの外部にあるか」を確認するだけで重なり判定ができます。
固定矩形Pに対するNFP(P,Q)の例。Qの参照点がNFP境界上にあるとき、PとQは接触（重ならない）状態になる。
矩形同士のNFP
矩形同士の場合、NFPは必ず矩形になります。PとQの左下座標・サイズをそれぞれ$(p_x, p_y, w_p, h_p)$、$(w_q, h_q)$とすると、NFP(P, Q)は以下の式で求まります。
\text{NFP}(P,Q)_\text{左下} = (p_x - w_q,\ p_y - h_q)
\text{NFP}(P,Q)_\text{サイズ} = (w_p + w_q) \times (h_p + h_q)
def calc_nfp(px, py, wp, hp, wq, hq):
"""固定矩形Pに対する移動矩形QのNFPを計算する。"""
return (px - wq, py - hq, wp + wq, hp + hq)
IFR（Inner-Fit Rectangle）
ビンからのはみ出しを防ぐためには、QをビンW×Hの内側に接するように動かしたときのQの参照点の軌跡、IFR（Inner-Fit Rectangle） を使います。
\text{IFR}_\text{左下} = (0,\ 0), \quad \text{IFR}_\text{サイズ} = (W - w_q) \times (H - h_q)
Qの参照点がIFRの範囲内にあれば、Qはビンからはみ出しません。
def calc_ifr(bin_w, bin_h_max, wq, hq):
"""ビン内でQの参照点が取れる範囲（IFR）を返す。"""
return (0.0, 0.0, bin_w - wq, bin_h_max - hq)
境界条件の重要性
実装上の重要なポイントとして、NFPの境界上は接触であり重なりではないという性質を正しく扱う必要があります。
y座標がNFPの内部にあるかを判定する箇所は、境界を含まない開区間で行う必要があります。
# 誤り: 境界上も「重なり」と判定してしまう
if ny <= y <= ny + nh:
# 正しい: 境界上を除外した開区間
if ny + 1e-9 < y < ny + nh - 1e-9:
同様にx方向の判定も開区間で行います。
# 誤り
if bx0 <= x <= bx1:
# 正しい
if bx0 + 1e-9 < x < bx1 - 1e-9:
この境界条件を誤ると、本来配置できる場所を「配置不可」と誤判定してしまい、正しいBL点が見つからなくなります。
Find2D-BL：走査線アルゴリズム
考え方
NFPを使うことで「Qの参照点がNFPの外部にあるか」という重なり判定が効率的にできるようになりました。これを利用してBL点を高速に探すアルゴリズムがFind2D-BLです。
X軸に平行な走査線をy=0から上方向に動かしながら、「その走査線上でNFPに覆われていない最小のx座標」を探します。最初に見つかった点がBL点です。
具体的には以下の手順で進みます。
- IFRを計算し、Qの参照点が取れるy範囲を得る
- 各既配置矩形のNFPを計算する
- y候補（IFR下端・各NFPの下端と上端）を昇順に並べる
- 各y候補において、IFRのx範囲内でNFPに覆われていない最小xを探す
- 最初に見つかった点（y最小 → x最小）がBL点
Pythonによる実装
def find_bl_point(wq, hq, placed_rects, bin_w, bin_h_max=1e9):
"""走査線アルゴリズムでBL点を探す。"""
ifr_x, ifr_y, ifr_w, ifr_h = calc_ifr(bin_w, bin_h_max, wq, hq)
if ifr_w < -1e-9 or ifr_h < -1e-9:
return None
# 各既配置矩形のNFPを計算
nfps = [
calc_nfp(px, py, pw, ph, wq, hq)
for (px, py, pw, ph) in placed_rects
]
# y候補を列挙してIFR範囲内に絞り込む
y_candidates = [ifr_y]
for (nx, ny, nw, nh) in nfps:
y_candidates.append(ny)
y_candidates.append(ny + nh)
y_upper = ifr_y + ifr_h
y_candidates = sorted(set(
y for y in y_candidates
if ifr_y - 1e-9 <= y <= y_upper + 1e-9
))
# 各y候補でNFPに覆われていない最小xを探す
for y in y_candidates:
active_x_intervals = [
(nx, nx + nw)
for (nx, ny, nw, nh) in nfps
if ny + 1e-9 < y < ny + nh - 1e-9
]
x = _find_min_x_not_covered(ifr_x, ifr_x + ifr_w, active_x_intervals)
if x is not None:
return (x, y)
return None
def _find_min_x_not_covered(x_min, x_max, blocked):
"""[x_min, x_max]内でblockedに含まれない最小xを返す。"""
x = x_min
improved = True
while improved:
improved = False
for (bx0, bx1) in blocked:
if bx0 + 1e-9 < x < bx1 - 1e-9:
x = bx1
improved = True
break
return x if x <= x_max + 1e-9 else None
計算量の改善
y候補はO(n)個あり、各y候補に対してNFPのx区間をスキャンするのにO(n)かかります。これだけ見るとO(n²)ですが、y候補を事前にソートしておくことでO(n log n)に改善できます。結果として全体の計算量はO(n²log n)になります。
| BL点探索 | 全体 | |
|---|---|---|
| 単純実装 | O(n³) | O(n⁴) |
| NFP版 | O(n log n) | O(n²log n) |
ベンチマーク結果
単純版とNFP版のBL法について、矩形数nを変えながら実行時間を計測しました。なお、両手法は全く同じ配置結果を出力することを確認しています。
実行時間の比較
| n | 単純版 (s) | NFP版 (s) | 高速化倍率 |
|---|---|---|---|
| 5 | 0.0001 | 0.0000 | 6x |
| 10 | 0.0008 | 0.0001 | 12x |
| 20 | 0.0100 | 0.0003 | 35x |
| 30 | 0.0483 | 0.0008 | 63x |
| 50 | 0.2928 | 0.0023 | 125x |
| 80 | 1.7243 | 0.0080 | 216x |
| 100 | 4.2966 | 0.0175 | 246x |
n=100において単純版が約4.3秒かかるのに対し、NFP版は約0.018秒と246倍の高速化を達成しています。
log-logプロット
log-logプロットで実測値と理論的な計算量の参照線を重ねると、単純版はO(n⁴)、NFP版はO(n²logn)の傾きに沿っていることが確認できます。
※ グラフはノートブック02_nfp_bottom_left.ipynbで再現できます。
配置結果の一致確認
高速化によって配置結果が変わっていないことを確認するため、ランダムに生成した50問×4ソート順（面積降順・幅降順・高さ降順・順不同）の計200ケースで両手法の出力を比較しました。全ケースで座標が完全に一致することを確認しています。
mismatch = 0
for seed in range(50):
rng = np.random.default_rng(seed=seed)
n = rng.integers(3, 20)
rects = [(float(rng.integers(2,8)), float(rng.integers(2,8))) for _ in range(n)]
bw = float(rng.integers(15, 30))
for key in ['area', 'width', 'height', 'none']:
ps, _ = bl_method(rects, bw, sort_key=key)
pn, _ = bl_method_nfp(rects, bw, sort_key=key)
for a, b in zip(ps, pn):
if abs(a[0]-b[0]) > 1e-6 or abs(a[1]-b[1]) > 1e-6:
mismatch += 1
print(f'不一致: {mismatch}件') # 0件
おわりに
本記事では、2次元矩形パッキング問題に対するBL法の実装と、NFPを用いた高速化について解説しました。
まとめ
- BL法は矩形を1つずつBL点（最も下で最も左）に配置する構築型ヒューリスティクスで、単純実装はO(n⁴)の計算量を持つ
- NFPを使うことで重なり判定を効率化でき、走査線アルゴリズム（Find2D-BL）と組み合わせることでO(n²logn)に改善できる
- 実装上の重要なポイントとして、NFPの境界上は接触であり重なりではないため、判定を開区間で行う必要がある
- n=100において単純版比246倍の高速化を達成しつつ、全く同じ配置結果を出力することを確認した
次回予告
次回は焼きなまし法を用いて配置順序を最適化し、充填率を改善する手法について解説する予定です。BL法単体では平均86.1%だった充填率が、焼きなまし法との組み合わせで平均94.9%（+8.8%）に改善できることを示します。
実装コードはGitHubで公開しており、StreamlitによるインタラクティブなデモURLも公開しています。ぜひ触ってみてください。