😄

自己紹介のグループ分けを数理最適化で解いてみた

2024/12/21に公開

この記事はJij Inc. Advent Calendar 2024の23日目の記事です。

はじめに

弊社Jij Inc.は基本的にリモートワークですが年に一度、ほぼ全メンバーがオフラインで集まる合宿というイベントを開催しています。2024年は11月6日から8日の3日間で開催されました。この合宿の初日に自己紹介というイベントがありました。今年入ったメンバーやリモートでしかやり取りをしたことがないメンバーも多いのでまずは自己紹介でアイスブレイクです。去年の合宿ではまだそれほど人数が多くなかったこともあり、全体で一人ずつ自己紹介していく方式でした。しかし今回の合宿では参加メンバーが29人と、一人ずつ自己紹介をしていくには人数が多いです。そのため、「5グループに分かれての自己紹介を4セッション繰り返す」という方式が取られました。この方式により少人数での自己紹介を通した密なやりとりが期待できます。一方で全メンバーの自己紹介を聞くことは出来なくなります。29人を5グループに分けるため1グループは多くて6人。セッション数は4なので最高でも20人のメンバーとしか同じグループになれません。

今回のグループ分けは担当者の方が人手でささっと決めたとのことでした。我々は数理最適化の企業なので社内行事にも積極的に数理最適化を使っていきたい!なので今回は数理最適化でより良いグループ分けを模索してみました。

結論

  • 非商用ソルバーSCIPで愚直に解くだけだと高速に高品質な解を得ることは出来なかった。
    • 5分の時間制限のもとでは実際のグループ分け(人手でささっと決めた)に劣る品質の解
  • しかし求解方法を工夫することで実際のグループ分けよりも優れたグループ分けを3分程度で得ることが出来た。
    • 工夫: 前のセッションから順番に最適化していく方法
  • 商用ソルバーGurobiを利用すると愚直に解くだけでも5分以内に実際のグループ分けよりも優れたグループ分けを得ることが出来た。
  • 商用ソルバーGurobiと、前のセッションから順番に最適化していく手法を組み合わせることで、さらに品質がいい解を得ることが出来た。
  • 注意点として、SCIPでもGurobiでも最適性の保証がある解は得られなかった。
    • ギャップが閉じず全然計算が終わらない...。なので上述の「5分程度で解が得られた」と言っているのは、5分の時間制限で求解を打ち切りその時点の解を取得しているということ

使用するライブラリ

Python環境に以下をインストールすることで本記事のコードを動かせます。ご興味持たれた方はぜひ試してみてください。

requirements.txt
jijmodeling == 1.10.1
ommx == 1.5.1
ommx-pyscipopt-adapter == 1.5.1
$ pip install -r requirements.txt

実際に使用されたグループ分けの確認

人手でささっと決めたという実際のグループ分けを確認し、どのような課題があるかを検討してみます。

実際のグループ分け
[
    [0, 5, 10, 15, 20, 25],  # グループA
    [1, 6, 11, 16, 21, 26],  # グループB
    [2, 7, 12, 17, 22, 27],  # グループC
    [3, 8, 13, 18, 23, 28],  # グループD
    [4, 9, 14, 19, 24],  # グループE
]  # セッション1

[
    [0, 9, 14, 18, 22, 28],  # グループA
    [1, 5, 10, 19, 23, 27],  # グループB
    [2, 6, 11, 15, 24, 25],  # グループC
    [3, 7, 12, 16, 20, 26],  # グループD
    [4, 8, 13, 17, 21],  # グループE
]  # セッション2

[
    [0, 8, 12, 16, 24, 27],  # グループA
    [1, 9, 13, 17, 20, 28],  # グループB
    [2, 5, 14, 18, 21, 26],  # グループC
    [3, 6, 10, 19, 22],  # グループD
    [4, 7, 11, 15, 23, 25],  # グループE
]  # セッション3

[
    [0, 7, 13, 19, 21, 26],  # グループA
    [1, 8, 14, 15, 22, 25],  # グループB
    [2, 9, 10, 16, 23, 28],  # グループC
    [3, 5, 11, 17, 24],  # グループD
    [4, 6, 12, 18, 20, 27],  # グループE
]  # セッション4

このように5グループに分かれての自己紹介を4回繰り返しました。各グループのリストに含まれる数字はメンバーのIDだと考えてください。私はこのグループ分けを分析してみて以下の2つを解決すべきと考えました。

  • 同じグループになれた他メンバーの人数について、メンバー間で差がある。
    • 20人の他メンバーと同じグループになることが出来たメンバーがいる。
      • 例えばメンバー0。メンバー0も含めて3人の該当メンバーがいる。
    • 一方で15人の他メンバーとしか同じグループになれなかったメンバーがいる。
      • 例えばメンバー11。メンバー11も含めて4人の該当メンバーがいる。
  • メンバー15とメンバー25は4セッション全てで同じグループになってしまっている。

この2点を解決するような数理モデルを検討していきます。

数理モデルの定式化

課題が解決されたグループ分けを得るべく、以下のような数理モデルを作成します。

  • 同じグループになることが出来た他メンバーの人数について、メンバーごとに出来るだけ差がないようにする。
    • 実際のグループ分けだと最も少ない人では15人の他メンバーとしか一緒のグループになることが出来なかった一方で最も多い人では20人の他メンバーと一緒のグループになることが出来ていて、5人の差があります。この差を縮めたいです。
    • このような状況では「最小値を最大化する」方法がよく使われます。つまり実際のグループ分けでいうと最小値の15人という数値を押し上げるというイメージです。
  • 全てのセッションで同じメンバーと同じグループという状況は禁止する。

定数・添字

N: メンバーの総数
M: グループ数
K: セッション数

i: メンバーを表す添字(0\leq i < N
j: グループを表す添字(0\leq j < M
k: セッションを表す添字(0\leq k < K

U: 同じ人と同じグループになっていい回数の上限(今回はU=K-1と設定します)
\text{group\_size}_{jk}: セッションkにおけるグループjの人数

決定変数

x_{ijk}: バイナリ変数。セッションkにおいてメンバーiがグループjに属するなら1、そうでないなら0

z: 同じグループになることが出来たメンバー数の最小値。もう少しちゃんと言うとメンバーiが同じグループになることが出来たメンバー数をz_iとするとz=\min\ z_i。実際のグループ割り当てで言えばz=15

x_{ijk}zだけでも要件通りの数理モデルを構築することは出来ますが、決定変数について二次の数理モデルになってしまいます(そのはずです)。数理モデルは基本的には線形に定式化した方がソルバーにとって解きやすいです。線形定式化のために追加で以下2つの決定変数を導入します。

y_{i_1i_2}: メンバーi_1とメンバーi_2が1度でも同じグループに属したら1、そうでないなら0

w_{i_1i_2jk}: セッションkにおいてメンバーi_1とメンバーi_2が同じグループjに属するなら1、そうでないなら0

目的関数

\max\ z

同じグループになることが出来たメンバー数の最小値を最大化します。
実際に使われたグループ分けではz=15です。数理最適化によってz=16以上、つまり1番少ないメンバーでも16人の他メンバーと同じグループになることが出来るようなグループ分けを得ることが出来れば本記事の取り組みは成功と言えるでしょう。

制約式

\sum_jx_{ijk}=1\ \ \forall\ i,k: セッションkにおいてメンバーiが属するグループは一つとする。

\sum_ix_{ijk}=\text{group\_size}_{jk}\ \ \forall\ j,k: セッションkでグループjに属する人数は\text{group\_size}_{jk}人とする。

\sum_{j,k}w_{i_1i_2jk}\leq U\ \ \forall i_1,i_2\ |\ i_1 < i_2: メンバーi_1とメンバーi_2が一緒になっていい回数はU回までとする。今回はU=K-1と設定する。つまり同じメンバーと全てのセッションで同じグループになることを禁止する制約となる。

x_{i_1jk}\geq w_{i_1i_2jk}\ \ \forall\ i_1,i_2,j,k\ |\ i_1<i_2: x_{ijk}w_{i_1i_2jk}の値を整合させるための制約。

x_{i_2jk}\geq w_{i_1i_2jk}\ \ \forall\ i_1,i_2,j,k\ |\ i_1<i_2: x_{ijk}w_{i_1i_2jk}の値を整合させるための制約。

x_{i_1jk}+x_{i_2jk}\leq w_{i_1i_2jk}+1\ \ \forall\ i_1,i_2,j,k\ |\ i_1<i_2: x_{ijk}w_{i_1i_2jk}の値を整合させるための制約。

\sum_{j,k}w_{i_1i_2jk}\geq y_{i_1,i_2}\ \ \forall\ i_1,i_2,j,k\ |\ i_1<i_2: w_{i_1i_2jk}y_{i_1i_2}の値を整合させるための制約。

z\leq\sum_{i_1,i_2:\ i_1=i\ \text{or}\ i_2=i,\ i_1<i_2}y_{i_1i_2}\ \ \forall\ i: y_{i_1i_2}zの値を整合させるための制約。

制約式についての補足

整合させるための制約がずらっと並び数理最適化に慣れていない方は少しびっくりしてしまうかもしれないので簡単に補足します。
決定変数x_{ijk}の定義は「セッションkにおいてメンバーiがグループjに属するなら1、そうでないなら0」でした。またw_{i_1i_2jk}の定義は「セッションkにおいてメンバーi_1とメンバーi_2がグループjに属するなら1、そうでないなら0」でした。これらの定義から、x_{ijk}w_{i_1i_2jk}の許される値の組み合わせが制限されます。

x_{i_1jk} x_{i_2jk} w_{i_1i_2jk} 定義上許されるか?
0 0 0
0 0 1 ×
0 1 0
0 1 1 ×
1 0 0
1 0 1 ×
1 1 0 ×
1 1 1

このテーブルの値を3つの「x_{ijk}w_{i_1i_2jk}の値を整合させるための制約」に代入してみてください。◯のケースは全て制約式を満たし、×のケースは全て制約式を満たさないことが確認出来るかと思います。

その他

実際のグループ分けを見て気付かれた方もいるかもしれませんが、全てのセッションにおいて各グループのメンバー1人はリーダーとして固定されています。数理モデルにも以下の制約を課します。

  • グループAにメンバー0は確定で加入
  • グループBにメンバー1は確定で加入
  • グループCにメンバー2は確定で加入
  • グループDにメンバー3は確定で加入
  • グループEにメンバー4は確定で加入

JijModelingによるモデリング

弊社のモデリングツールJijModelingで、こちらの数理モデルを書いたものが以下です。

import jijmodeling as jm

N = jm.Placeholder("N", description="number of members")
M = jm.Placeholder("M", description="number of groups")
K = jm.Placeholder("K", description="number of sessions")
U = jm.Placeholder("U", description="max number of same person in same group")

group_size_ph = jm.Placeholder("group\_size", ndim=2, description="group size")
fix_leader_ph = jm.Placeholder("fix\_leader", ndim=1, description="fix leader in each group")

x = jm.BinaryVar("x", shape=(N, M, K), description="セッションkにおいてメンバーiがグループjに属するなら1、そうでないなら0")
y = jm.BinaryVar("y", shape=(N, N), description="メンバーi\_1とメンバーi\_2が1度でも同じグループに属したら1、そうでないなら0")
w = jm.BinaryVar("w", shape=(N, N, M, K), description="セッションkにおいてメンバーi\_1とメンバーi\_2がグループjに属するなら1、そうでないなら0")
z = jm.ContinuousVar("z", description="同じグループになることが出来たメンバー数の最小値", lower_bound=0, upper_bound=N-1)

i1 = jm.Element("i1", belong_to=(0, N))
i2 = jm.Element("i2", belong_to=(0, N))
i3 = jm.Element("i3", belong_to=(0, N))
j = jm.Element("j", belong_to=(0, M))
k = jm.Element("k", belong_to=(0, K))

problem = jm.Problem("group_self_introduction", sense=jm.ProblemSense.MAXIMIZE)

problem += z

# あるセッションで一人のメンバーはただ一つのグループに所属
problem += jm.Constraint("assign one group", jm.sum(j, x[i1, j, k]) == 1, forall=[i1, k])

# グループサイズの指定
problem += jm.Constraint("group size", jm.sum(i1, x[i1, j, k]) == group_size_ph[k, j], forall=[j, k])

# 同じ人と同じグループになれる回数の制限
problem += jm.Constraint("max grouping with same person", jm.sum([j, k], w[i1, i2, j, k]) <= U, forall=[i1, (i2, i1 < i2)])

# x、w、y、zの整合性
problem += jm.Constraint("assign same group 1", x[i1, j, k] >= w[i1, i2, j, k], forall=[i1, (i2, i1 < i2), j, k])
problem += jm.Constraint("assign same group 2", x[i2, j, k] >= w[i1, i2, j, k], forall=[i1, (i2, i1 < i2), j, k])
problem += jm.Constraint("assign same group 3", x[i1, j, k] + x[i2, j, k] <= w[i1, i2, j, k] + 1, forall=[i1, (i2, i1 < i2), j, k])
problem += jm.Constraint("assign same group 4", jm.sum([j, k], w[i1, i2, j, k]) >= y[i1, i2], forall=[i1, (i2, i1 < i2)])
problem += jm.Constraint("z is min of other member", z <= jm.sum([i2, (i3, (i2 < i3) & ((i2==i1) | (i3==i1)))], y[i2, i3]), forall=i1)

# リーダーは固定
problem += jm.Constraint("fix leader", x[fix_leader_ph[j], j, k] == 1, forall=[j, k])

problem


出力
JijModelingでは数式に近い書き方でモデリングを行うことができ、また画像のように定義した数理モデルをLaTeX表示することで、自身が作成した数理モデルをチェックすることが出来ます。

以下のようにインスタンスデータも定義しておきます。JijModelingのコードの、対応するPlaceholderに代入されます。

instance_data = {'N': 29,
 'M': 5,
 'K': 4,
 'U': 3,
 'group\\_size': [[6, 6, 6, 6, 5],
  [6, 6, 6, 5, 6],
  [6, 6, 5, 6, 6],
  [6, 5, 6, 6, 6]],
 'fix\\_leader': [0, 1, 2, 3, 4],
}

求解

数理モデルが出来たので求解していきます。

SCIPによる求解

フリーの数理最適化ソルバーとして定評があるSCIPソルバーで今回の問題を解いてみます。

愚直な求解

以下のコードでSCIPを使った求解が出来ます。特に工夫や設定などを行わずに解いています。

from ommx_pyscipopt_adapter import instance_to_model, model_to_solution

# problem、instance_dataからOMMXのインスタンスを生成
interpreter = jm.Interpreter(instance_data)
ommx_instance = interpreter.eval_problem(problem)

# OMMXのインスタンスからSCIPのモデルを生成
model = instance_to_model(ommx_instance)

# ログを出力する
model.hideOutput(False)

# SCIPで最適化
model.optimize()

# SCIPの解をOMMXのSolutionに変換
ommx_solution_by_scip = model_to_solution(model, ommx_instance)

if ommx_solution_by_scip.feasible:
    print("制約を満たす解が得られました")
    print(f'目的関数値: {ommx_solution_by_scip.objective}')
else:
    print("制約を満たす解が得られませんでした")

このコードを実行してみると...計算が全然終わりません。今回の問題は後段で出てくる商用ソルバーGurobiに10時間以上解かせても計算が終わらないので、SCIPやGurobiなどのソルバーでパッと厳密な最適解が得られる問題ではないようです。

工夫1: 求解時間を制限する

SCIPやGurobiなどの厳密ソルバーを使うときに求解時間に制限を設けるのは簡単かつ重要な工夫です。厳密ソルバーは今得られている解が最適であるという保証が得られるまで計算を続けます。そのため計算が全く終わらない時でも、実はかなり優れた解が得られている場合があります。

今回は求解時間制限を300秒としましょう。300秒が経過したら、最適性の保証が得られていなくてもソルバーの計算処理が終了します。300秒とした理由は「今回くらいの規模のグループ分けで5分以上待ちたくないなー」という私の感覚によります。

先ほどのコードにmodel.setParam("limits/time", 300)を追加すればOKです。

from ommx_pyscipopt_adapter import instance_to_model, model_to_state

# problem、instance_dataからOMMXのインスタンスを生成
interpreter = jm.Interpreter(instance_data)
ommx_instance = interpreter.eval_problem(problem)

# OMMXのインスタンスからSCIPのモデルを生成
model = instance_to_model(ommx_instance)

# 時間制限を設定
model.setParam("limits/time", 300)

# ログを出力する
model.hideOutput(False)

# SCIPで最適化
model.optimize()

# SCIPの解をOMMXのSolutionに変換
ommx_state_by_scip = model_to_state(model, ommx_instance)
ommx_solution_by_scip = ommx_instance.evaluate(ommx_state_by_scip)

if ommx_solution_by_scip.feasible:
    print("制約を満たす解が得られました")
    print(f'目的関数値: {ommx_solution_by_scip.objective}')
else:
    print("制約を満たす解が得られませんでした")
# ログは省略
目的関数値: 13.0

目的関数値、すなわちz=13...。これは4セッションを通して一緒のグループになれた他メンバーの人数が13人のメンバーがいるということを意味します。実際のグループ分けでは最小でも15人だったので、これは実際のグループ分けに劣る結果だということになります。

工夫2: 前のセッションから順番に解く

次の工夫は「前のセッションから順番に解く」というものです。具体的には以下です。

  1. 1セッション目だけの最適化問題を解く
  2. 2セッション目までの最適化問題を解く。ただし1セッション目の決定変数は1の解で固定する
  3. 3セッション目までの最適化問題を解く。ただし1,2セッション目の決定変数は1,2の解で固定する
  4. 4セッション目までの最適化問題(つまり全体)を解く。ただし1,2,3セッション目の決定変数は1,2,3の解で固定する

この方法により1回の最適化計算の負担を減らし、全体としてより品質の良い解が得られる可能性があります。ただし最終的に得られた解が元の問題の最適解である保証が無いという点には注意が必要です。こちらの手法でも1-4全体でソルバーの求解時間制限を300秒までとして求解してみます。

以下にコードを載せます。数理モデルにも修正が入っている(前の解の数値を固定するための制約、fix_x・fix_y・fix_wの部分)ため改めて全体を載せておきます。

import jijmodeling as jm
import numpy as np
import time
from ommx_pyscipopt_adapter import instance_to_model, model_to_solution

# 各種設定
num_members = 29
num_groups = 5
num_sessions = 4
time_limit_sec = 300

N = jm.Placeholder("N", description="number of members")
M = jm.Placeholder("M", description="number of groups")
K = jm.Placeholder("K", description="number of sessions")
U = jm.Placeholder("U", description="max number of same person in same group")

group_size_ph = jm.Placeholder("group\_size", ndim=2, description="group size")
fix_leader_ph = jm.Placeholder("fix\_leader", ndim=1, description="fix leader in each group")

x = jm.BinaryVar("x", shape=(N, M, K), description="セッションkにおいてメンバーiがグループjに属するなら1、そうでないなら0")
y = jm.BinaryVar("y", shape=(N, N), description="メンバーi\_1とメンバーi\_2が1度でも同じグループに属したら1、そうでないなら0")
w = jm.BinaryVar("w", shape=(N, N, M, K), description="セッションkにおいてメンバーi\_1とメンバーi\_2がグループjに属するなら1、そうでないなら0")
z = jm.ContinuousVar("z", description="同じグループになることが出来たメンバー数の最小値", lower_bound=0, upper_bound=N-1)

i1 = jm.Element("i1", belong_to=(0, N))
i2 = jm.Element("i2", belong_to=(0, N))
i3 = jm.Element("i3", belong_to=(0, N))
j = jm.Element("j", belong_to=(0, M))
k = jm.Element("k", belong_to=(0, K))

problem = jm.Problem("group_self_introduction", sense=jm.ProblemSense.MAXIMIZE)

problem += z

# あるセッションで一人のメンバーはただ一つのグループに所属
problem += jm.Constraint("assign one group", jm.sum(j, x[i1, j, k]) == 1, forall=[i1, k])

# グループサイズの指定
problem += jm.Constraint("group size", jm.sum(i1, x[i1, j, k]) == group_size_ph[k, j], forall=[j, k])

# 同じ人と同じグループになれる回数の制限
problem += jm.Constraint("max grouping with same person", jm.sum([j, k], w[i1, i2, j, k]) <= U, forall=[i1, (i2, i1 < i2)])

# x、w、y、zの整合性
problem += jm.Constraint("assign same group 1", x[i1, j, k] >= w[i1, i2, j, k], forall=[i1, (i2, i1 < i2), j, k])
problem += jm.Constraint("assign same group 2", x[i2, j, k] >= w[i1, i2, j, k], forall=[i1, (i2, i1 < i2), j, k])
problem += jm.Constraint("assign same group 3", x[i1, j, k] + x[i2, j, k] <= w[i1, i2, j, k] + 1, forall=[i1, (i2, i1 < i2), j, k])
problem += jm.Constraint("assign same group 4", jm.sum([j, k], w[i1, i2, j, k]) >= y[i1, i2], forall=[i1, (i2, i1 < i2)])
problem += jm.Constraint("z is min of other member", z <= jm.sum([i2, (i3, (i2 < i3) & ((i2==i1) | (i3==i1)))], y[i2, i3]), forall=i1)

# リーダーは固定
problem += jm.Constraint("fix leader", x[fix_leader_ph[j], j, k] == 1, forall=[j, k])

# 変数の値を固定する。前のセッションから順番に解く方法で利用する
fix_x_ph = jm.Placeholder("fix\_x", ndim=3, description="fix x in some session")
fix_y_ph = jm.Placeholder("fix\_y", ndim=2, description="fix y in some session")
fix_w_ph = jm.Placeholder("fix\_w", ndim=4, description="fix w in some session")
problem += jm.Constraint("fix x", x[i1, j, k] == fix_x_ph[i1, j, k], forall=[i1, j, (k, fix_x_ph[i1, j, k] > -0.99)])
problem += jm.Constraint("fix y", y[i1, i2] == fix_y_ph[i1, i2], forall=[i1, (i2, (fix_y_ph[i1, i2] > -0.99) & (i1 < i2))])
problem += jm.Constraint("fix w", w[i1, i2, j, k] == fix_w_ph[i1, i2, j, k], forall=[i1, (i2, i1 < i2), j, (k, fix_w_ph[i1, i2, j, k] > -0.99)])


# 最初のiter(1セッション目だけを最適化)ではfixなし(-1が入っている場合は固定を行わない)
fix_x = np.full((num_members, num_groups, num_sessions), -1)
fix_y = np.full((num_members, num_members), -1)
fix_w = np.full((num_members, num_members, num_groups, num_sessions), -1)

# 合計求解時間
total_solver_elapsed_time = 0

# 残りのイテレーション数
rest_iter = num_sessions

# 前から1セッションずつ最適化を行う。次のセッションの求解をする時、それより前のセッションはそれまでの解に固定する
for iter_count in range(num_sessions):
    # インスタンスデータの生成
    num_sessions_iter = iter_count + 1
    num_sessions_nextiter = iter_count + 2
    instance_data = {'N': num_members,
    'M': num_groups,
    'K': num_sessions_iter,
    'U': num_sessions - 1,
    'group\\_size': [[6, 6, 6, 6, 5],
    [6, 6, 6, 5, 6],
    [6, 6, 5, 6, 6],
    [6, 5, 6, 6, 6]],
    'fix\\_leader': [0, 1, 2, 3, 4],
    'fix\\_x': fix_x,
    'fix\\_y': fix_y,
    'fix\\_w': fix_w
    }

    # problem、instance_dataからOMMXのインスタンスを生成
    interpreter = jm.Interpreter(instance_data)
    ommx_instance = interpreter.eval_problem(problem)

    # OMMXのインスタンスからSCIPのモデルを生成
    model = instance_to_model(ommx_instance)

    # 今回の求解に使える時間を計算する
    time_limit_sec_iter = (time_limit_sec - total_solver_elapsed_time) / rest_iter

    # 今回の求解時間を計測する
    solver_start_time = time.perf_counter()

    # 時間制限を設定
    model.setParam("limits/time", time_limit_sec_iter)

    # ログを出力する
    model.hideOutput(False)

    # SCIPで最適化
    model.optimize()

    # 合計求解時間を更新
    total_solver_elapsed_time += time.perf_counter() - solver_start_time

    # SCIPの解をOMMXのSolutionに変換
    ommx_solution_by_scip = model_to_solution(model, ommx_instance)

    # 次のiteration向けのfix_x, fix_y, fix_wを生成(既に得られたセッションは固定する)
    fix_x = np.full((num_members, num_groups, num_sessions), -1)
    fix_y = np.full((num_members, num_members), -1)
    fix_w = np.full((num_members, num_members, num_groups, num_sessions), -1)
    for _, dv in ommx_solution_by_scip.decision_variables.iterrows():
        if dv["name"].item() == "x":
            i, j, k = dv["subscripts"].item()
            if dv["value"].item() > 0.99:
                fix_x[i, j, k] = 1
            else:
                fix_x[i, j, k] = 0
        elif dv["name"].item() == "y":
            i1, i2 = dv["subscripts"].item()
            if dv["value"].item() > 0.99:
                fix_y[i1, i2] = 1
        elif dv["name"].item() == "w":
            i1, i2, j, k = dv["subscripts"].item()
            if dv["value"].item() > 0.99:
                fix_w[i1, i2, j, k] = 1
            else:
                fix_w[i1, i2, j, k] = 0
    
    rest_iter -= 1

# ログは省略
目的関数値: 16.0

目的関数値が16の解を得ることが出来ました!これは4セッションを通して一緒のグループになれた他メンバーの人数が、1番少ない人でも16人であるようなグループ割り当てを得ることが出来たということです。これは実際のグループ割り当ての15人よりも優れた数値です。ちなみに求解時間は201.8秒でした。

SCIPでそのまま数理モデルを解くと、5分の時間制限では人手で作った実際のグループ割り当てを上回る解を得ることが出来ませんでしたが、「前のセッションから順番に解く」という工夫をすることで、実際のグループ割り当てを超える解を3分程度で得ることが出来ました。このように一度にソルバーに投げる問題サイズを小さくする工夫は実際の数理最適化プロジェクトにおいても有効な手法です。

Gurobiによる求解

せっかくなので商用ソルバーのGurobiでも解いてみました。コードは省略させて頂きます。結果は実際に使われたもの、SCIPの解も併せて、以下のテーブルにまとめました。

手法 所要時間 [sec] 目的関数値(一緒になれたメンバー数の最小値) [人]
人手(実際のグループ分け) - 15
SCIPナイーブ 300 13
SCIP前から法 201.8 16
Gurobiナイーブ 300 16
Gurobi前から法 41.5 17
  • ナイーブ → 問題全体をそのままソルバーで求解
  • 前から法 → 「前のセッションから順番に解く」方法で求解

双方、全体でソルバーの求解時間は300秒以内の制限を設けています。

Gurobiはナイーブな方法でも人手のグループ分けを上回る品質の解が得られており、さらに「前から法」と組み合わせることでより優れた解を高速に得ることが出来ています。
一つ注意点として、SCIPとGurobiで動かしたマシンが異なるため公平なベンチマークではないという点に言及しておきます。

  • SCIPよりもGurobiの方が高スペックなマシンで動かしている
    • SCIPは私の普段使いのマシンで動かしている
    • Gurobiのライセンスは高額でありインストールできるマシン数にも制限がある。そのため性能を引き出すためにハイスペックなVMにインストールされている。

SCIPもGurobiのマシンで動かして公平なベンチマークを取ろうか...と迷いましたが、「気軽にやるとこんなもん」という結果も欲しいなと思って今回の形を取りました(両方試せよと言われたらそれはそう)。

Gurobi前から法で得られたグループ分けを詳細に確認してみると、以下のようになっていることが確認出来ました(確認用のコードなどは省略させて頂きます)。

同じグループになることが出来たメンバー数の最小値(目的関数値): 17
同じグループになることが出来たメンバー数の最大値: 18

これは得られたグループ分けにおいて、各メンバーが他メンバーと同じグループになることが出来た人数は17人、もしくは18人であるということです。数理最適化によってメンバー間の不平等が少ない非常に良いグループ分けを得ることが出来ました。

最後に

本記事では自己紹介のグループ分けに数理最適化を適用し、人手で得たグループ分けよりも優れた解を得ることが出来ました(少なくとも今回設定した目的関数においては)。

本編では触れませんでしたが、当然SCIP・Gurobiで得られた解はいずれも「同じメンバーと全てのセッションで同じグループになることは禁止する」という制約を満たしています。グループ分けの担当者が今回の解を確認した時に他にも満たしてほしい条件があれば、それも制約式として追加していけるのは数理最適化の良いところです。

個人的にブラッシュアップしたい点としては、各メンバーが他メンバーに対して持つ「話してみたい度」「自己紹介で交流してみたい度」みたいなものを考慮することです。例えば普段同じチームに属して作業しているメンバーや、年次が高いメンバー同士はお互い既によく知っていて、自己紹介で同じグループになることの効用が小さいかもしれません。そのような条件も反映したらより良いグループ分けが出来るかもしれません。

募集

Jijでは数理最適化・量子コンピュータに関連するソフトウェア開発を行っています。
OSSも複数開発、提供しています。Discordのコミュニティへぜひ参加してください!
https://discord.gg/Km5dKF9JjG
また,Jijでは各ポジションを絶賛採用中です!
ご興味ございましたらカジュアル面談からでもぜひ!ご連絡お待ちしております。
数理最適化エンジニア(採用情報
量子コンピューティングソフトウェアエンジニア(採用情報
Rustエンジニア(アルゴリズムエンジニア)(採用情報
研究チームPMO(採用情報
オープンポジション等その他の採用情報はこちら

Discussion