生産現場で役に立つかもしれない最適化手法4【GLPK】(工数の平準化2)

[<--目次へ](https://qiita.com/ki073/items/508b3c20eb03078703dd)

## はじめに
仕様や原理などは後日加筆予定です。

工数の平準化するプログラムを作ってみました。出荷予定を指定することで生産計画を作ることきます。プログラムは２つに分かれていまので、順番に実行してください。
GLPKのインストールや使い方は上の目次リンクからたどってください。
ファイルは全てUTF-8で作ってください。LinuxやMacOSではファイル名の制限はありませんが、Windowsでは下記の通り日本語を含むファイル名は使えませんので注意してください。

###Windowsのgusekを使うときの注意(2021.8.3 追記)

プログラムとデータファイルは１つの同じフォルダに入れてください。
プログラムの中で読み書きするファイルの名前に日本語が含まれる場合は
csv_driver: unable to open ???.csv - No such file or directory
のようにファイルがみつからないとのエラーになります。必ず英数字にしてください。フィルダ名には日本語を使っても大丈夫なようです。
文字化けがある時には、上の目次リンクからたどって「インストールと簡単な使い方」を見てください。WindowsのExcelで作られるcsvファイルはShift-JISで出力されるようです。Excel-2016からUTF-8出力を選べるような記事を見たのですが、手元にあったExcel-2016では無理でした。Windows10のメモ帳でShift-JISからUTF-8に変換できます。

## プログラム１(2021.8.3 読み書きするファイルの名称を修正)
出荷計画ファイル(shipping_schedule.csv)が必要ですので、下記のようなCSV形式のファイルを表計算ソフトなどで作成しておいてください。
製品在庫上限(max_stock_quantity)と製品在庫下限(min_stock_quantity)はデータブロックに書き込むことで製品単位で指定可能です。
プログラム中に
param ~ default 0;
のように書かれている場合はデータブロックで変更可能です。データブロックに何も書かなければdefault値が使われます。

```GLPK:manufacuring31.mod
# version 0.1.1
# 出荷計画ファイルから読み込むデータ
set Product_Date dimen 2; # CSVファイルから読み込む製品名と出荷日
set Product := setof{(p, d) in Product_Date}p; # 製品名を取り出す
set Date := {0..9999} inter setof{(p, d) in Product_Date}d; # 工場稼働日を取り出し、昇順に並び替える
# 日付の範囲が0..9999を超える時は変更のこと
param shipping_quantity{Product, Date} default 0; # 出荷数
# ファイルを読み込む
table table_shipping IN "CSV" "shipping_schedule.csv":
	Product_Date <- [Product, Date], 
	shipping_quantity ~ Shipping;

## 以下２つのパラメータはデータブロックで変更可能 ##
# 製品ごとの在庫下限を設定
param min_stock_quantity{Product} default 0;
# 製品ごとの在庫上限を設定 出荷量の最大値の２倍としている
param max_stock_quantity{p in Product} default 2*max{d in Date}shipping_quantity[p,d];

# 稼働日ごとの出荷後の在庫下限を求める
param min_stock_quantity2{p in Product, d in Date}:= min_stock_quantity[p];
# 稼働日ごとの出荷後の在庫上限を求める
param max_stock_quantity2{p in Product, d in Date}:= max_stock_quantity[p];

# ファイルから読み込んだ情報を表示
printf("\n出荷計画\n");
printf{d in Date} "\t%s", d; printf("\t合計\n");
for{p in Product}{
	printf "%s", p;
	printf{d in Date} "\t%g", shipping_quantity[p, d];
	printf "\t%g\n", sum{d in Date}shipping_quantity[p, d];
}
printf("\n");
printf("在庫設定\n");
printf("\t下限\t上限\n");
for{p in Product}{
	printf "%s", p;
	printf "\t%g\t%g\n", min_stock_quantity[p], max_stock_quantity[p];
	
}
printf("\n\n");

# プログラム2で読み込む中間ファイルを書き出す
table table_result{p in Product, d in Date} OUT "CSV" "shpping_stock_inter.csv":
	p ~ Product, d ~ Date,  shipping_quantity[p, d] ~ Shipping,
	min_stock_quantity2[p, d] ~ MinStock, max_stock_quantity2[p, d] ~ MaxStock;

data;
# 製品在庫上限を変更
param max_stock_quantity :=
	"製品5" 100;
# 製品在庫下限を変更
param min_stock_quantity :=
	"製品4" 3
	"製品2" 10;
end;

```

以下の出荷計画ファイル(shipping_schedule.csv)が必要ですので作成しておいてください。製品名、出荷日、出荷数の順です。Excelでも作成可能です。

このファイルから、製品名と工場稼働日データを読み込みます。出荷日として１度も現れない日は休業日として判断されます。どこかに出荷数0として入れておいてください。(４日がないので製品4,4,0を追加)

出荷計画はないが在庫を増やしたい場合は、製品名と期間内の適当な日付、出荷数0と書き込んでおいてください。ただし、出荷数０の場合は製品在庫上限をが最大出荷数の２倍、すなわち０となりますので後でエラーが出ます。データブロックで製品在庫上限を在庫を増やす数（以上）を指定しておいてください。(製品5)

```shipping_schedule.csv
Product,Date,Shipping
製品1,1,100
製品1,2,100
製品1,3,200
製品1,5,100
製品1,8,100
製品1,9,100
製品1,10,100
製品1,11,100
製品1,12,100
製品2,1,200
製品2,2,200
製品2,3,200
製品2,5,400
製品2,8,200
製品2,9,200
製品2,10,200
製品2,11,200
製品2,12,200
製品3,2,300
製品3,6,300
製品3,8,300
製品3,10,300
製品3,12,300
製品4,2,2
製品4,3,2
製品4,5,2
製品4,9,2
製品4,11,2
製品4,4,0
製品5,1,0

```

##プログラム2 (2021.8.3 読み書きするファイルの名称を修正)
上のプログラムで作成されたshpping_stock_inter.csvを読み込んで、生産計画などを出力します。
日ごとの最大工数(max_man_hour)、日ごとの最低工数(min_man_hour)は必ずデータブロックに書きてください。消さないように。
それらが適切でないと、"PROBLEM HAS NO PRIMAL FEASIBLE SOLUTION"と表示され、結果は出力されません。
「日平均工数」が出力されますのでそれを参考に書き換えてください。

1バッチの製造数(batch_quantity):
バッチ単位の製造を想定していますので、標準的なサイズがあればプログラム中のdefault値(10としています)を書き換えて、個別に変更したい時にはデータブロックで指定してください。
 
1単位(個)を製造するための工数(man_hour_for_manufacturing):
default値を1としていますので、１個作るのに1工数かかるとしています。製品ごとに書き換えがデータブロックでできます。

```GLPK:manufacuring32.mod
# version 0.1.1
# 出荷在庫計画ファイルから読み込むデータ
set Product_Date dimen 2; # 読み込む製品と日付
set Product := setof{(p, d) in Product_Date}(p); # 製品名を取り出す
set Date := setof{(p, d) in Product_Date}(d); # 工場稼働日を取り出す
param min_stock_quantity0{Product, Date} default 0; # 読み込む在庫下限
param max_stock_quantity{Product, Date} default 0; # 読み込む在庫上限
param shipping_quantity{Product, Date} default 0;  # 読み込む出荷数
# ファイル読み込み
table table_shipping IN "CSV" "shpping_stock_inter.csv" :
	Product_Date <- [Product, Date], 
	shipping_quantity ~ Shipping, min_stock_quantity0 ~ MinStock, max_stock_quantity ~ MaxStock;

# 在庫計画ファイルから読み込むデータ
set ProductInStockfile; # 読み込む製品名
param initial_stock2{ProductInStockfile} default 0; # 読み込む初期在庫
param final_stock2{ProductInStockfile} default 0; # 読み込む最終在庫
table table_stock IN "CSV" "stock_plan.csv" :
	ProductInStockfile <- [Product], 
	initial_stock2 ~ InitialStock,
	final_stock2 ~ FinalStock;
# 出荷在庫計画ファイルに現れる製品名の情報を取り出す。在庫計画ファイルになければ0をする
param initial_stock{p in Product} := if (p in ProductInStockfile) then initial_stock2[p] else 0;
param final_stock{p in Product} := if (p in ProductInStockfile) then final_stock2[p] else 0;

# 計算で使用する日付情報を求める
set DatesForCumulate{d in Date} := {w in Date : d>=w}; # 累積計算用日付Set
set FirstDay := setof{d in Date : card(DatesForCumulate[d])=1}(d); # 初日
set LastDay := setof{d in Date : card(DatesForCumulate[d])=card(Date)}(d); # 最終日

# 初期在庫設定が在庫上限を超える時はエラーにする
check {p in Product, d in FirstDay}: initial_stock[p]<=shipping_quantity[p,d]+max_stock_quantity[p,d];


## 以下２つのパラメータはデータブロックで必ず設定すること ##
param max_man_hour{Date}; # 日ごとの最大工数
param min_man_hour{Date}; # 日ごとの最低工数
## 次の２つのパラメータはデータブロックで変更可能 ##
param batch_quantity{Product} default 10; # 1バッチの製造数
param man_hour_for_manufacturing{Product} default 1; # 1単位(個)を製造するための工数

param subtotal_of_shipping_quantity{p in Product} := sum{d in Date}shipping_quantity[p, d]; # 製品別出荷数の合計
param stock_increase{p in Product} := final_stock[p]-initial_stock[p]; # 在庫増

# 最終日の在庫下限を設定した最終日の在庫数に書き換え
param min_stock_quantity{p in Product, d in Date} :=
	if {d} within LastDay then final_stock[p] else min_stock_quantity0[p, d];

# ファイルから読み込んだ情報を表示
printf "出荷計画\n";
printf{d in Date} "\t%s", d; printf("\t在庫増\t合計\n");
for{p in Product}{
	printf "%s", p;
	printf{d in Date} "\t%g", shipping_quantity[p,d];
	printf "\t%g\t%g\n", stock_increase[p], stock_increase[p]+sum{d in Date}shipping_quantity[p,d];
}
printf "工数計";
printf{d in Date}: "\t%g", sum{p in Product}shipping_quantity[p,d]*man_hour_for_manufacturing[p];
printf "\t%g", sum{p in Product}stock_increase[p]*man_hour_for_manufacturing[p];
printf "\t%g\n", sum{p in Product}((sum{d in Date}shipping_quantity[p,d])+stock_increase[p])*man_hour_for_manufacturing[p];
printf "最終在庫は在庫下限であるため"; printf "製造数等は表示より増加することがあります\n";
printf "日平均工数\t%g\n", (sum{p in Product}((sum{d in Date}shipping_quantity[p,d])+stock_increase[p])*man_hour_for_manufacturing[p])/card(Date);
printf "\n";
printf "在庫設定\n";
printf("\t初期\t最終\n");
printf{p in Product} "%s\t%g\t%g\n", p, initial_stock[p], final_stock[p];
printf "\n";

# 累積出荷数
param cumulated_shipping_quantity{p in Product, d in Date} := sum{d1 in DatesForCumulate[d]}shipping_quantity[p, d1];
# 累積製造数下限上限
param stock_lower_limit_include_cumulated_shipping{p in Product, d in Date} := cumulated_shipping_quantity[p,d]+min_stock_quantity[p,d]-initial_stock[p];
param stock_upper_limit_include_cumulated_shipping{p in Product, d in Date} := cumulated_shipping_quantity[p,d]+max_stock_quantity[p,d]-initial_stock[p];

# 累積製造数下限上限をバッチ数に変換
param stock_lower_limit_batches{p in Product, d in Date} integer := max(ceil(stock_lower_limit_include_cumulated_shipping[p, d]/batch_quantity[p]), 0);
param stock_upper_limit_batches{p in Product, d in Date} integer := max(ceil(max(stock_upper_limit_include_cumulated_shipping[p, d], stock_lower_limit_include_cumulated_shipping[p, d])/batch_quantity[p]), 0);

# 製品、稼働日へのバッチ数を割り当てする変数
var assign_n_batches{Product, Date} integer >=0;

# 各稼働日の工数上限の制約
s.t. upper_limit_man_hours{d in Date}: sum{p in Product}assign_n_batches[p, d]*batch_quantity[p]*man_hour_for_manufacturing[p]<=max_man_hour[d];
# 各稼働日の工数下限の制約
s.t. lower_limit_man_hours{d in Date}: sum{p in Product}assign_n_batches[p, d]*batch_quantity[p]*man_hour_for_manufacturing[p]>=min_man_hour[d];
# 各稼働日の累積製造数+在庫数の上限で制約
s.t. upper_batch_limit{p in Product, d in Date}: sum{d1 in DatesForCumulate[d]}assign_n_batches[p,d1]<=stock_upper_limit_batches[p,d];
# 各稼働日の累積製造数+在庫数の下限で制約
s.t. lower_batch_limit{p in Product, d in Date}: sum{d1 in DatesForCumulate[d]}assign_n_batches[p,d1]>=stock_lower_limit_batches[p,d];
# 累積製造数+在庫数の下限との差を最小化
minimize minimize_stock: sum{p in Product, d in Date}(assign_n_batches[p,d]-stock_lower_limit_batches[p,d]);

solve;

# 製造数
param manufactured_quantity{p in Product, d in Date} := assign_n_batches[p,d]*batch_quantity[p];
# 工数
param manufactured_man_hour{p in Product, d in Date} := manufactured_quantity[p,d]*man_hour_for_manufacturing[p];
# 累積出荷数
param cumulated_shipped_quantity{p in Product, d in Date} := sum{d1 in DatesForCumulate[d]}shipping_quantity[p,d1];
# 初期在庫+累積製造数
param cumulated_manufactured_quantity{p in Product, d in Date} := initial_stock[p]+sum{d1 in DatesForCumulate[d]}manufactured_quantity[p,d1];

# 結果を出力
printf "\n製造計画\n";
printf{d in Date} "\t%s", d; printf("\t合計\n");
for{p in Product}{
	printf "%s", p;
	printf{d in Date} "\t%g", manufactured_quantity[p,d];
	printf "\t%g\n", sum{d in Date}manufactured_quantity[p,d];
}
printf "工数計";
printf{d in Date}: "\t%g", sum{p in Product}manufactured_man_hour[p,d];
printf "\t%g\n", sum{p in Product, d in Date}manufactured_man_hour[p,d];

# 在庫見込み
printf "\n在庫見込み(出荷日は出荷後)\n";
printf "\t初期";
printf{d in Date} "\t%s", d; printf("\t最終日設定\n");
for{p in Product}{
	printf "%s", p;
	printf "\t%g", initial_stock[p];
	printf{d in Date} "\t%g", cumulated_manufactured_quantity[p,d]-cumulated_shipped_quantity[p,d];
	printf "\t%g\n", final_stock[p];
}

data;
# 日ごとの最大工数 必須
param max_man_hour default 500
	4 450; # 4日だけ特別に指定
# 日ごとの最低工数 必須
param min_man_hour default 400;
# 1バッチの数量を変更
param batch_quantity :=
	"製品4" 5
	"製品2" 100
	"製品1" 50;
end;

```

製品の期間前の在庫と期間後の在庫(こちらは最低在庫となっているのでそれより在庫が増える可能性があります)を設定した在庫計画ファイル(stock_plan.csv)が必要です。
出荷計画にない製品は書かれていても無視され、抜けている製品は在庫0として扱われます。

```stock_plan.csv
Product,InitialStock,FinalStock
製品1,100,200
製品2,100,300
製品4,2,3
製品5,0,100

```

###動作確認
GLPK 5.0, MacOS 10.14.6
GLPK 5.0, CentOS8
gusek 0-2-24(64bits), Windows10