数理最適化でアドベントカレンダーのお菓子を決めよう

去年のアドベントカレンダーの最適化の記事を参考に問題を作成しました。

（記事のデータや条件を少し変えています）

問題

16種類、136個のお菓子があります。
このお菓子をアドベントカレンダー用の25個の袋に入れてください。

データはCSVに入っています。

advent.csv

商品,小分け数,単位,大きさ
A. ブルボン プチチョコラングドシャ,1,1,90
B. グリコ 冬のくちどけポッキー,2,1,90
C. グリコ 冬のきらめきポッキー,2,1,90
D. 亀田製菓 亀田のうす焼えび,3,1,90
E. キタノ商事 トップオブザポップ バター,1,1,90
F. 江崎グリコ プリッツ(ロースト 塩バター),2,1,90
G. ロッテ トッポ(ザ・ショコラ),2,1,90
H. グリコ 生チーズのチーザ カマンベールチーズ,1,1,90
I. 明治 メルティーキッスプレミアムショコラ,14,2,5
J. 明治 メルティーキッスフルーティー濃いちご,14,2,5
K. 明治 メルティーキッス初摘み濃抹茶,14,2,5
L. 亀田製菓 ハッピーターン,24,2,10
M. ネスレ日本 キットカット,12,2,10
N. 亀田製菓 亀田のまがりせんべい,10,2,15
O. 国産小麦の厚切りバウムクーヘン,10,2,15
P. ミックスゼリー(135g),24,2,5

商品	小分け数	単位	大きさ
A. ブルボン プチチョコラングドシャ	1	1	90
B. グリコ 冬のくちどけポッキー	2	1	90
C. グリコ 冬のきらめきポッキー	2	1	90
D. 亀田製菓 亀田のうす焼えび	3	1	90
E. キタノ商事 トップオブザポップ バター	1	1	90
F. 江崎グリコ プリッツ(ロースト 塩バター)	2	1	90
G. ロッテ トッポ(ザ・ショコラ)	2	1	90
H. グリコ 生チーズのチーザ カマンベールチーズ	1	1	90
I. 明治 メルティーキッスプレミアムショコラ	14	2	5
J. 明治 メルティーキッスフルーティー濃いちご	14	2	5
K. 明治 メルティーキッス初摘み濃抹茶	14	2	5
L. 亀田製菓 ハッピーターン	24	2	10
M. ネスレ日本 キットカット	12	2	10
N. 亀田製菓 亀田のまがりせんべい	10	2	15
O. 国産小麦の厚切りバウムクーヘン	10	2	15
P. ミックスゼリー(135g)	24	2	5

条件その1

お菓子は大きいもの（列単位が1）と小さいもの（列単位が2）があります。
小さいものは2人で分けられるように2の倍数の個数にしてください。
小さいものの小分け数は偶数になっています。
また、いろいろ食べられるように1種類は1袋に4個までとします。

条件その2

袋の大きさは90なので、列大きさの合計が90までしか入りません。
また、大きさの合計は90 * 25なので、25個とも90になります。

条件その3

前の日と同じお菓子が入っていると残念な気持ちになるので、同じ商品が2日連続しないようにしてください。

解答

このような問題は数理最適化で解くことができます。
数理最適化では、数理モデルを作成してソルバーで解を求めます。

詳しくは下記を参考にしてください。

まずは、数理モデルを考えます。以降ではお菓子を商品とも呼ぶことにします。

数理モデル

数理モデルは、変数、目的関数、制約条件から構成されます。

変数

変数は、2種類考えます。

1つは、あるお菓子をある日に選ぶかどうかを表す変数です。これは0-1変数を使います（V01とします）。
もう1つは、何単位買うかを表す変数です。これは整数変数を使います（_Viとします）。

それぞれ16種類 x 25日分用意します。

ここで、_Viに単位を掛けたものが個数になります。これをVarとします。Varは数式ですが、便宜上、これも変数とみなしましょう。

目的関数

今回は、目的関数はありません。

制約条件

3つの制約条件は、次のように考えます。

• V01を1にすればVarを4個まで可能
• 各日の大きさの和は90
• 同じ商品を2日連続にしない

また、商品は小分け数までしかないので、下記の制約条件も必要です。

• 商品ごとの小計は小分け数

まとめると次のようになります。

構成要素	内容
変数	V01: 入れるかどうか（商品ごと、日ごと） _Vi: 単位数（商品ごと、日ごと） Var: 個数 （商品ごと、日ごと）
目的関数	なし
制約条件	V01を1にすればVarを4個まで可能 各日の大きさの和は90 同じ商品を2日連続にしない 商品ごとの小計は小分け数

実装

PolarsとPython−MIPを使って数理モデルを作成して解いてみましょう。

Polarsを使ったPython-MIPの数理モデルの作り方は次を参照してください。

import numpy as np
import polars as pl
from mip import Model, xsum

def xmul(s1, s2):
    """内積"""
    return np.array(s1) * np.array(s2)

df = pl.read_csv("advent.csv").join(
    pl.DataFrame(range(1, 26), ["日"]), how="cross"
)
# 数理モデル
m = Model()
V01 = m.add_var_tensor((len(df),), "V01", var_type="B")
_Vi = m.add_var_tensor((len(df),), "_Vi", var_type="I")
# 変数追加
df = df.with_columns(V01=V01, Var=xmul(_Vi, df["単位"]))

for row in df.iter_rows(named=True):
    # V01を1にすればVarを最大個数まで可能
    m += row["Var"] <= 4 * row["V01"]

for _, group in df.group_by("日", maintain_order=True):
    # 各日の大きさの和は90
    m += xsum(xmul(group["Var"], group["大きさ"])) == 90

for _, group in df.group_by("商品", maintain_order=True):
    for i in range(1, 25):
        # 同じ商品を2日連続にしない
        m += xsum(group.filter(pl.col("日").is_between(i, i + 1))["V01"]) <= 1

for _, group in df.group_by("商品", maintain_order=True):
    # 商品ごとの小計は小分け数
    m += xsum(group["Var"]) == group["小分け数"][0]

m.verbose = 0
# ソルバー実行
m.optimize()
print(m.status)  # 結果のステータス
# 結果の列を追加
df = df.with_columns(
    Val=df["Var"].to_numpy().astype(float).round().astype(int)
)
# 結果の表示
print(
    df.filter(df["Val"] > 0)
    .group_by("日")
    .agg(pl.struct("商品", "Val"))
    .sort("日")
)

出力

OptimizationStatus.OPTIMAL

shape: (25, 2)

日	商品
i64	list[struct[2]]
1	[{"I. 明治 メルティーキッスプレミアムショコラ",4}, {"K. 明治 メルティーキッス初摘み濃抹茶",2}, {"L. 亀田製菓 ハッピーターン",4}, {"P. ミックスゼリー(135g)",4}]
2	[{"J. 明治 メルティーキッスフルーティー濃いちご",4}, {"M. ネスレ日本 キットカット",4}, {"N. 亀田製菓 亀田のまがりせんべい",2}]
3	[{"I. 明治 メルティーキッスプレミアムショコラ",4}, {"K. 明治 メルティーキッス初摘み濃抹茶",2}, {"L. 亀田製菓 ハッピーターン",4}, {"P. ミックスゼリー(135g)",4}]
4	[{"G. ロッテ トッポ(ザ・ショコラ)",1}]
5	[{"K. 明治 メルティーキッス初摘み濃抹茶",4}, {"L. 亀田製菓 ハッピーターン",2}, {"O. 国産小麦の厚切りバウムクーヘン",2}, {"P. ミックスゼリー(135g)",4}]
6	[{"J. 明治 メルティーキッスフルーティー濃いちご",4}, {"M. ネスレ日本 キットカット",4}, {"N. 亀田製菓 亀田のまがりせんべい",2}]
7	[{"K. 明治 メルティーキッス初摘み濃抹茶",4}, {"L. 亀田製菓 ハッピーターン",2}, {"O. 国産小麦の厚切りバウムクーヘン",2}, {"P. ミックスゼリー(135g)",4}]
8	[{"B. グリコ 冬のくちどけポッキー",1}]
9	[{"I. 明治 メルティーキッスプレミアムショコラ",2}, {"N. 亀田製菓 亀田のまがりせんべい",2}, {"O. 国産小麦の厚切りバウムクーヘン",2}, {"P. ミックスゼリー(135g)",4}]
10	[{"D. 亀田製菓 亀田のうす焼えび",1}]
11	[{"B. グリコ 冬のくちどけポッキー",1}]
12	[{"F. 江崎グリコ プリッツ(ロースト 塩バター)",1}]
13	[{"H. グリコ 生チーズのチーザ カマンベールチーズ",1}]
14	[{"D. 亀田製菓 亀田のうす焼えび",1}]
15	[{"I. 明治 メルティーキッスプレミアムショコラ",4}, {"J. 明治 メルティーキッスフルーティー濃いちご",4}, {"M. ネスレ日本 キットカット",2}, {"O. 国産小麦の厚切りバウムクーヘン",2}]
16	[{"G. ロッテ トッポ(ザ・ショコラ)",1}]
17	[{"F. 江崎グリコ プリッツ(ロースト 塩バター)",1}]
18	[{"C. グリコ 冬のきらめきポッキー",1}]
19	[{"L. 亀田製菓 ハッピーターン",4}, {"M. ネスレ日本 キットカット",2}, {"O. 国産小麦の厚切りバウムクーヘン",2}]
20	[{"D. 亀田製菓 亀田のうす焼えび",1}]
21	[{"E. キタノ商事 トップオブザポップ バター",1}]
22	[{"A. ブルボン プチチョコラングドシャ",1}]
23	[{"J. 明治 メルティーキッスフルーティー濃いちご",2}, {"K. 明治 メルティーキッス初摘み濃抹茶",2}, {"L. 亀田製菓 ハッピーターン",4}, {"N. 亀田製菓 亀田のまがりせんべい",2}]
24	[{"C. グリコ 冬のきらめきポッキー",1}]
25	[{"L. 亀田製菓 ハッピーターン",4}, {"N. 亀田製菓 亀田のまがりせんべい",2}, {"P. ミックスゼリー(135g)",4}]

全条件を満たした解が得られました。

以上