CPLEXサンプル ProdPlan（モデル１生産台数と部品調達先の計画）を読み解く


CPLEXのサンプルのProdPlan（生産台数と部品調達先の計画）を読み解きます。

<サンプル導入ディレクトリ>opl\examples\opl\ProdPlan\totalprod.mod
<サンプル導入ディレクトリ>opl\examples\opl\ProdPlan\period.mod
にあります。

* 実行環境
    * CPLEX 22.1.1
    * Windows 11 64bit

# 問題

コンピューター組み立て工場の3期分の生産計画を作成します。この工場は、MultimediaBusiness、MultimediaHome、および BasicHome という 3 つのコンピューター・モデルを組み立てることができます。
それぞれのモデル毎、期ごとに最小需要数と最大需要数、そして在庫可能数があります。
各期ごとに組立可能な台数があります。
各モデル毎にCPU やモデムなどの必要なコンポーネントは異なります。各コンポーネントは内製、外注それぞれでコストが異なります。

モデルを2段階に分けて生産計画を作っていきます。この記事はモデル２を対象にします。

モデル１：組み立てる各タイプのコンピューターの総数と、必要なコンポーネントを内製、外注する数を決定します。

[モデル２：各期間に組み立てる数、販売する数、および在庫数を決定します。](https://qiita.com/kawada2017/items/d3f4b7279bf03c053d47)


![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/223763/45636047-1539-05cc-2da6-ca25c0bc6ce6.png)

* サンプル

●コード

https://github.com/hkwd/230823cplexsample/blob/master/231124ProdPlan/totalprod_hkwd.mod
https://github.com/hkwd/230823cplexsample/blob/master/231124ProdPlan/period_hkwd.mod

この記事でのサンプルは、オリジナルのサンプルとほぼ同じ内容ですが、日本語コメントを入れたり順序を入れ替えたり、添え字を入れ替えたりして読みやすくしています

●サンプルデータ

https://github.com/hkwd/230823cplexsample/blob/master/231124ProdPlan/production_hkwd.dat

https://github.com/hkwd/230823cplexsample/blob/master/231124ProdPlan/period_hkwd.dat



## 利用データ
1.  各期毎の最大組み立て可能数

|期間|値|
|:----|:----|
|W1|35|
|W2|20|
|W3|25|


２.  コンピュータータイプのコンポーネントや需要量

|コンピュータタイプ|コンポ―ネント|価格|最大可能在庫数|最小需要W1|最小需要W2|最小需要W3|最大需要W1|最大需要W2|最大需要W3|
|:----|:----|:----|:----|:----|:----|:----|:----|:----|:----|
|MultimediaBusiness|{CPU250 Modem56K SoundCard VideoCard CDRW}|8000|5|5|6|7|15|16|17|
|MultimediaHome|{CPU250 Modem56K SoundCard VideoCard}|4000|5|4|5|6|14|15|16|
|BasicHome|{CPU150 Modem288 VideoCard}|2000|5|3|4|5|13|14|15|

３.各コンポーネントのサプライヤーと価格スロープ
発注数が21を超えると割引単価になります。

|コンポーネント|サプライヤー|20個までの単価|21個以上の単価|
|:----|:----|:----|:----|
|CPU150|CheapComponents|500|450|
|CPU150|DiscoutCPUs|475|450|
|CPU250|CheapComponents|325|300|
|CPU250|DiscoutCPUs|400|250|
|Modem288|CheapComponents|125|100|
|Modem288|ModemsRUs|150|75|
|Modem56K|CheapComponents|90|90|
|Modem56K|ModemsRUs|95|85|
|SoundCard|CheapComponents|50|50|
|SoundCard|InHouse|55|45|
|VideoCard|CheapComponents|75|75|
|VideoCard|InHouse|80|50|
|CDRW|CheapComponents|105|90|
|CDRW|CDsForASong|100|100|


## 決定変数
0. 各コンピューターを全部で何台ずつ作るか
0. 必要な各コンポーネント数
0. どのコンポーネントをどのサプライヤーから何個仕入れるか


## 目的関数
粗利益＝コンピューターの売上-コンポーネントのコスト


## 制約
0.  各コンピューター数は最大可能組み立て数以下
0.  各コンピューター数は最大需要数以下
0.  各コンピューター数は最小需要数以上
0.  組み立てるコンピュータータイプの台数と必要なコンポーネント数が一致する
0.  供給されるコンポーネント数と必要なコンポーネント数が一致する

## 問題の種類
混合整数計画

# 解
先に解を紹介します。

## 目的関数
粗利益の合計:442,340

## 決定変数
１. 各コンピューターを全部で何台ずつ作るか

|コンピュータタイプ|値|価格|売上|
|:----|:----|:----|:----|
|MultimediaBusiness|48|8000|384,000|
|MultimediaHome|20|4000|80,000|
|BasicHome|12|2000|24,000|
| | | |488,000|

２. 必要な各コンポーネント数
３. どのコンポーネントをどのサプライヤーから何個仕入れるか

|コンポ―ネント|サプライヤー|発注数|20個までの単価|21個以上の単価|20個までのコスト|21個以上のコスト|総コスト|
|:----|:----|:----|:----|:----|:----|:----|:----|
|CPU150|CheapComponents|0|500|450|0|0|0|
|CPU150|DiscoutCPUs|12|475|450|5700|0|5,700|
|CPU250|CheapComponents|0|325|300|0|0|0|
|CPU250|DiscoutCPUs|68|400|250|8000|12000|20,000|
|Modem288|CheapComponents|12|125|100|1500|0|1,500|
|Modem288|ModemsRUs|0|150|75|0|0|0|
|Modem56K|CheapComponents|0|90|90|0|0|0|
|Modem56K|ModemsRUs|68|95|85|1900|4080|5,980|
|SoundCard|CheapComponents|0|50|50|0|0|0|
|SoundCard|InHouse|68|55|45|1100|2160|3,260|
|VideoCard|CheapComponents|0|75|75|0|0|0|
|VideoCard|InHouse|80|80|50|1600|3000|4,600|
|CDRW|CheapComponents|48|105|90|2100|2520|4,620|
|CDRW|CDsForASong|0|100|100|0|0|0|
| | | | | | | |45,660|




# OPLコード解説
OPLコードの解説を行っていきます。

## 利用データ

まずコンピュータータイプとコンポーネントのデータの定義を行っています。
データの中身の定義はproduction_hkwd.datファイルの中で行っています。

```java:コンピュータータイプとコンポーネントの定義(mod)
{string} ComputerTypes = ...;
{string} ComponentTypes = ...;
```


```java:コンピュータータイプとコンポーネントの定義(dat)
ComputerTypes = { MultimediaBusiness, MultimediaHome, BasicHome };
ComponentTypes = { CPU150, CPU250, Modem288, Modem56K, SoundCard, VideoCard, CDRW };
```


各コンピュータタイプには、コンポーネントのセット、販売価格、および在庫可能最大数があります。

```java:コンピュータ(mod)
 *   --------------------------------------------------- */
tuple computersToBuild
{
  {string} components;
  int      price;
  int      maxInventory;
}

computersToBuild Computers[ComputerTypes] = ...;
```

```java:コンピュータ(dat)
Computers = [<{CPU250, Modem56K,SoundCard,VideoCard,CDRW},8000,5>,
                     <{CPU250,Modem56K,SoundCard,VideoCard},4000,5>,
                     <{CPU150,Modem288,VideoCard},2000,5>];

```

問題ブラウザで見ると以下のようになっています。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/223763/f5af3cc6-9d35-37dd-62d2-6d3ce5d5c838.png)

期間を定義しています。3週間と定義しています。

```java:期間(mod)
期間数
int NbPeriods = ...;
range Periods = 1..NbPeriods;
```
```java:期間(dat)
NbPeriods = 3;
```
問題ブラウザで見ると以下のようになっています。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/223763/eda21fa3-751a-53cf-a289-b1aa011c9046.png)

各コンピューターの種類には、各期間に販売できる最大数量と最小数量があります。 これらの値は、各コンピューターの種類の総組立台数を計算するために使用されます。 各期間の最大組立可能台数もあります。

```java:需要数、組立可能台数(mod)
//期ごとの最大の組立可能台数
int MaxBuildPerPeriod[1..NbPeriods] = ...;

//コンピュータの種類と期ごとの最小最大の需要量
int MinDemand[ComputerTypes][1..NbPeriods] = ...;
int MaxDemand[ComputerTypes][1..NbPeriods] = ...;
```
```java:需要数、組立可能数(dat)
MaxBuildPerPeriod = [35 20 25];

MinDemand = [[5,6,7],
         [4,5,6],
         [3,4,5]];

MaxDemand = [[15,16,17],
         [14,15,16],
           [13,14,15]];
```

問題ブラウザで見ると以下のようになっています。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/223763/125340b1-cbc8-111f-52d3-d398535b5774.png)
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/223763/d15d483a-267f-5809-20d1-fc3ba255ebef.png)
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/223763/fecd445c-810b-a56c-6c44-f6400b0fb8d1.png)

次に、[execute](https://www.ibm.com/docs/ja/icos/20.1.0?topic=keywords-execute)でOPLのスクリプトをつかってfor文で全期分の需要数および組立可能数を計算しています。このようにdatから読まずにデータ定義することも可能です。

```java:全期分の需要数、組立可能数(mod)
//全期分の最大組立数
int MaxBuild = 0;
//コンピュータの種類ごと全期分の最小最大の需要量
int TotalMaxDemand[ComputerTypes];
int TotalMinDemand[ComputerTypes];
//未使用float TotalBuild[ComputerTypes];
//for文で初期化
execute INITIALIZE {
  for(var p in Periods) {
    MaxBuild = MaxBuild + MaxBuildPerPeriod[p];
    for(var c in ComputerTypes) {
      TotalMaxDemand[c] = MaxDemand[c][p] + TotalMaxDemand[c];
      TotalMinDemand[c] = MinDemand[c][p] + TotalMinDemand[c];
    }
  }
}
```

問題ブラウザで見ると以下のようになっています。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/223763/b1e2b16b-1f5e-6fab-7ab0-1f3a1458dec1.png)
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/223763/dc3c1326-42e2-efdb-47d8-4b08868a0da2.png)
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/223763/54b5b626-8dc6-1e0d-df48-a6761c701ec8.png)

次に、部品サプライヤーの定義です。
```java:部品サプライヤー(mod)
{string} Suppliers = ...;
```
```java部品サプライヤー(dat)
Suppliers = { CheapComponents ModemsRUs DiscoutCPUs CDsForASong InHouse };
```

次に、部品サプライヤーのコストの定義です。
NbIncrementsやSupplyCostIncrementsは、20個までの発注と21個以上の発注でコストが1回変化するという定義を行っています。21個以上発注すると値引きがあることを意味しています。
ComponentsはComponentTypes毎にサプライヤーとコストの変化を定義しています。
ComponentSupplierMatchesはComponentsと同じ情報ですが、componentInfoのタプルの中身を展開しています。目的関数や制約を書きやすくするために別のセットを作っています。

```java:部品サプライヤーのコスト(mod)
//サプライヤーのコスト変化ポイントの定義
int NbCostChanges = ...;         // 区分線形のブレイク数　Number of "breaks" in PLF
int SupplyCostCostChanges[1..NbCostChanges] = ...;

//コンポーネントのタイプ毎、サプライヤーごとのコストの傾斜
tuple componentInfo
{
  string supplier;
  int    costSlope[1..NbCostChanges+1];   //CostChangesppと同じであるべき。should be CostChangespp
}
{componentInfo} Components[ComponentTypes] = ...;

tuple supplierMatch {
  string         component;
  componentInfo  componentInformation;
}
{supplierMatch} ComponentSupplierMatches = { <i,j> | i in ComponentTypes, j in Components[i] };
```
```java:部品サプライヤー(dat)
NbCostChanges = 1;
SupplyCostCostChanges = [20];
          
Components = [ {<CheapComponents, [500,450]>, <DiscoutCPUs, [475,450]>},
          {<CheapComponents, [325,300]>, <DiscoutCPUs, [400,250]>},
          {<CheapComponents, [125,100]>, <ModemsRUs, [150,75]>},
          {<CheapComponents, [90,90]>, <ModemsRUs, [95,85]>},
          {<CheapComponents,  [50,50]>, <InHouse, [55,45]>},
          {<CheapComponents,  [75,75]>, <InHouse, [80,50]>},
          {<CheapComponents,  [105,90]>, <CDsForASong, [100,100]>}];

```


問題ブラウザで見ると以下のようになっています。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/223763/5014ef8b-d006-5be6-9592-d99ce90682fd.png)
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/223763/5252e085-b242-90d5-964e-d77cbd370553.png)

## 決定変数
次に決定変数の定義を行っています。
dvarが決定変数を意味します。

```java:決定変数
//組立数 -- 各タイプのコンピューターを何台組み立てるか。
dvar int+ Build[ComputerTypes];
//必要な各コンポーネント数 -- Buildに基づく 
dvar int+ NecessaryComponents[ComponentTypes];
//コンポーネント供給数 -- 各コンポーネント タイプがどのサプライヤーによっていくつ供給されるか。
dvar int+ SuppliedComponents[ComponentSupplierMatches];
```

これらの決定変数は、求められた解では以下のようになりました。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/223763/fdde8d64-9ccb-10ee-0711-523b1f04226a.png)
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/223763/d91fc22a-8bc9-0dc2-7062-90367c4db011.png)
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/223763/fd2c2edd-3c1b-7ca2-62b2-79f93327534b.png)




## 目的関数
次に目的関数の定義を行っています。

関数の定義に先立って、決定式のコスト（Cost）、売上（Sales）、粗利益(GrossProfit)を作ります。

ます、コストを計算します。

[sumキーワード](https://www.ibm.com/docs/ja/icos/20.1.0?topic=keywords-sum)は関連式の集合の合計値を計算するキーワードです。
[piecewiseキーワード](https://www.ibm.com/docs/ja/icos/20.1.0?topic=keywords-piecewise)は区分線形関数を表現するための OPL キーワードです。
piecewiseキーワードをつかってサプライヤーのコストが発注数が20を超えると割引価格になることを定義しています。このようにある点で線形関数の傾きが変わる関数を区分線形関数といいます。


```java:コスト
dexpr float Cost =       
  sum(m in ComponentSupplierMatches)
    piecewise(i in 1..NbCostChanges) {
      m.componentInformation.costSlope[i] -> SupplyCostCostChanges[i];   
      m.componentInformation.costSlope[NbCostChanges+1]
    } SuppliedComponents[m];
```

発注数（SuppliedComponents[m]）が20個（SupplyCostCostChanges[1]）までは
m.componentInformation.costSlope[1]×SuppliedComponents[m]
で、21個からは
m.componentInformation.costSlope[2]×(SuppliedComponents[m]-20)
になるという定義になっています。

例えばCPU250をDiscoutCPUsというサプライヤーに発注した場合、costSlopeは20個までは400ドル、21個からは250ドルになります。ここに68個発注した場合、400ドル×20個+250ドル×(68-20)個=20,000ドルになるということです。

|コンポ―ネント|サプライヤー|発注数|20個までの単価|21個以上の単価|20個までのコスト|21個以上のコスト|総コスト|
|:----|:----|:----|:----|:----|:----|:----|:----|
|CPU150|CheapComponents|0|500|450|0|0|0|
|CPU150|DiscoutCPUs|12|475|450|5700|0|5,700|
|CPU250|CheapComponents|0|325|300|0|0|0|
|CPU250|DiscoutCPUs|68|400|250|8000|12000|20,000|
|Modem288|CheapComponents|12|125|100|1500|0|1,500|
|Modem288|ModemsRUs|0|150|75|0|0|0|
|Modem56K|CheapComponents|0|90|90|0|0|0|
|Modem56K|ModemsRUs|68|95|85|1900|4080|5,980|
|SoundCard|CheapComponents|0|50|50|0|0|0|
|SoundCard|InHouse|68|55|45|1100|2160|3,260|
|VideoCard|CheapComponents|0|75|75|0|0|0|
|VideoCard|InHouse|80|80|50|1600|3000|4,600|
|CDRW|CheapComponents|48|105|90|2100|2520|4,620|
|CDRW|CDsForASong|0|100|100|0|0|0|
| | | | | | | |45,660|

ここでは価格は一回しか変化しない、つまり傾きは一回しか変化しないという定義になっていますが、SupplyCostCostChanges[i]は配列になっており、複数回変化するように定義することもできます。例えば41個目からはさらに安くなるという定義も可能です。


次は売上です。各コンピューターの価格と組み立て数の掛け算です。

```java:売上
//売上 -- ビルドの線形関数
dexpr float Sales=
  sum(p in ComputerTypes) Computers[p].price * Build[p];
```
|コンピュータタイプ|値|価格|売上|
|:----|:----|:----|----:|
|MultimediaBusiness|48|8,000|384,000|
|MultimediaHome|20|4,000|80,000|
|BasicHome|12|2,000|24,000|
| | | |488,000|

次は目的変数となる粗利益です。売上-コストです。
これを最大化maximizeします。

```java:目的関数　粗利益
//粗利益
dexpr float GrossProfit = Sales - Cost;

////目的関数 粗利益の最大化
maximize GrossProfit;
```
## 制約
次に制約の定義を行っています。
subject toの中カッコで囲みます。
```python:制約（mod）
subject to{

---略----
}

```

### 制約1: 各コンピューターは最大組み立て数以下

ctXXX:で制約にラベルを付けています。

```java:制約1: 各コンピューターは最大組み立て数以下
ctPlantCapacity:
  sum(p in ComputerTypes) Build[p] <= MaxBuild;
```
問題ブラウザで見ると以下のようになります。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/223763/b72ffb22-d2f8-eaa7-2f93-7ecd4e252b0d.png)

図にすると以下のようなイメージになります。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/223763/9915d73b-3209-661b-8793-483b9e2ad040.png)

### 制約2、3:各コンピューターは最大需要数以下、最小需要数以上

forall(p in ComputerTypes)でコンピュータタイプ数分の制約を作っています。



```java:制約2、3:各コンピューターは最大需要数以下、最小需要数以上
forall(p in ComputerTypes){
  ctComputerTypeMaxDemand: Build[p] <= TotalMaxDemand[p];
  ctComputerTypeMinDemand: Build[p] >= TotalMinDemand[p];      
}  
```
問題ブラウザで見ると以下のようになります。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/223763/3f7bfa13-424f-0a9d-2d82-104a79471317.png)

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/223763/1bd806f0-4f80-fadc-96dc-65a2368dd1e5.png)


図にすると以下のようなイメージになります。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/223763/0bad33d2-9213-cd0f-c6ed-34b74d2bca09.png)
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/223763/08dc9b4f-10fe-d361-dc63-131e4160bda5.png)



### 制約4:組み立てるコンピュータータイプの台数と必要なコンポーネント数が一致する
forall(c in ComponentTypes)でコンポ―ネントタイプ数分の制約を作っています。
共通のコンポーネントを持つ（p in ComputerTypes: c in Computers[p].components）のBuild数の合計とNecessaryComponentsが一致することという制約です。

```java:制約4:組み立てるコンピュータータイプの台数と必要なコンポーネント数が一致する
forall(c in ComponentTypes)
  ctDetermineAmtNecessary:      
    sum(p in ComputerTypes: c in Computers[p].components) 
      Build[p] == NecessaryComponents[c];
```
問題ブラウザで見ると以下のようになります。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/223763/86c2c736-a8b5-9a74-c43e-d2c4e99b7f4b.png)



図にすると以下のようなイメージになります。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/223763/6b413b49-5e4e-b337-24f6-cd67f966cc7c.png)

### 制約5:供給されるコンポーネント数と必要なコンポーネント数が一致する
forall(c in ComponentTypes)でコンポ―ネントタイプ数分の制約を作っています。
各サプライヤーのコンポーネントを持つComputerTypes（m in ComponentSupplierMatches: c == m.component）のSuppliedComponents数の合計とNecessaryComponentsの一致することという制約です。

```java:制約5:供給されるコンポーネント数と必要なコンポーネント数が一致する
forall(c in ComponentTypes)
  ctDetermineAmtSupplied:      
    sum(m in ComponentSupplierMatches: c == m.component) 
       SuppliedComponents[m] == NecessaryComponents[c];

```
問題ブラウザで見ると以下のようになります。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/223763/d50d67d1-2b48-1bba-5e4a-1191cf4908ee.png)

図にすると以下のようなイメージになります。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/223763/00afa4ab-0de7-f8a2-0c91-2c5b7414610d.png)

## モデルの実行

ここまででモデルとしては完成していますが、モデルの実行は二つの実行をするためにOPLスクリプトでモデルを実行しています。

### モデルの実行１：生産台数と部品調達先の計画モデル
[mainキーワード](https://www.ibm.com/docs/ja/icos/20.1.0?topic=steps-defining-main-block)はフロー制御のスクリプト・ブロックを導入するための OPL キーワードです。このブロックにOPLスクリプトを記載します。
[thisOplModel](https://www.ibm.com/docs/ja/icos/20.1.0?topic=steps-loading-necessary-structures)でここまで書いてきたモデルを呼び出し、[generate()](https://www.ibm.com/docs/ja/icos/20.1.0?topic=steps-generating-optimization-model)で生成して、[cplex.solve()](https://www.ibm.com/docs/ja/icos/20.1.0?topic=steps-solving-current-optimization-model)で求解します。


```java:一つ目のモデルの求解
main {
   //数量決定モデルの設定
   var quantity = thisOplModel;
   quantity.generate();
   if (cplex.solve()) {
```

目的関数の結果をcplex.getObjValue()で得て出力しています。

```java:粗利益出力(mod)
writeln("Gross profit: ", cplex.getObjValue());
```

```java:粗利益出力（スクリプトログ）
Gross profit: 442340
```

quantity.Build[c].solutionValueでBuild決定変数にアクセスをしています。


```java:各コンピュータータイプ毎の組立台数の合計を算出して出力(mod)
for(c in quantity.ComputerTypes) {
  totalbuild[c] = quantity.Build[c].solutionValue;
  writeln("   ", c, ":  ", totalbuild[c]);
```

```java:各コンピュータータイプ毎の組立台数の合計を算出して出力（スクリプトログ）
   MultimediaBusiness:  48
   MultimediaHome:  20
   BasicHome:  12
```

内製で生産するコンポーネントを出力します。quantity.SuppliedComponents.componentInformation.supplier == "InHouse"でInHouseから供給されるコンポーネントを絞り込んでいます。

```java:内製で生産するコンポーネント(mod)

 writeln("Components to build in-house:");
 for(m in quantity.SuppliedComponents){
   if(quantity.SuppliedComponents[m] > 0 
   && m.componentInformation.supplier == "InHouse") {
     writeln("  ", m.component, ": ", quantity.SuppliedComponents[m].solutionValue);
}
```

```java:内製で生産するコンポーネント（スクリプトログ）
Components to build in-house:
  SoundCard: 68
  VideoCard: 80
```

### モデル２の実行：各期間最適化モデル

次にモデル２の各期間最適化モデルを解きます。period_hkwd.modとそのデータのperiod_hkwd.datを呼び出します。
このモデル２については[別記事](https://qiita.com/kawada2017/items/d3f4b7279bf03c053d47)で解説します。
main内で、[modファイルやdatファイルを読み出すことができます](https://www.ibm.com/docs/ja/icos/20.1.0?topic=steps-loading-necessary-structures)。

```java:期間最適化のモデル読込(mod)
var source = new IloOplModelSource("period_hkwd.mod");
var def = new IloOplModelDefinition(source);
var newCplex = new IloCplex();
var dates = new IloOplModel(def, newCplex);
var data = new IloOplDataSource("period_hkwd.dat");
dates.addDataSource(data);
data = dates.dataElements;
```

この時点では、datファイルのTotalBuildには以下のダミーデータ入っています。
```java:モデル１の最適化結果でデータを書き換え前(period_hkwd.dat)
//一つ目の問題の解を入力するのでダミー値
TotalBuild = [0 0 0];
```

これをモデル１で得たtotalbuild[c]、つまり以下の値で書き換えます。
   MultimediaBusiness:  48
   MultimediaHome:  20
   BasicHome:  12

```java:モデル１の最適化結果でデータを書き換える。(mod)
for(c in dates.ComputerTypes) {
  data.TotalBuild[c] = totalbuild[c];
}
```

モデル２の求解をして結果を書き出します。

各期事各コンピュータータイプ毎の組立台数、販売台数、在庫台数が出力されます。
なお、モデル２の結果や定義は問題ブラウザで見ることはできません。
このようなスクリプトで結果を出力する必要があります。

```java:モデル２の求解をして結果を書き出し(mod)
dates.generate();
if (newCplex.solve()) {
  // 出力　Output
  for(var p in dates.Periods) {
    writeln("Period: ", p)
    for(c in dates.ComputerTypes) {
      writeln("  ", c, ":  ")
      writeln("      Build: ", dates.Build[c][p].solutionValue)
      writeln("      Sell:  ", dates.Sell[c][p].solutionValue)
      writeln("      Hold:  ", dates.InStockAtEndOfPeriod[c][p].solutionValue)
    }
  }
}
dates.end();
```

```java:モデル２の求解をして結果を書き出し(スクリプトログ)
Period: 1
  MultimediaBusiness:  
      Build: 20
      Sell:  15
      Hold:  5
  MultimediaHome:  
      Build: 7
      Sell:  4
      Hold:  3
  BasicHome:  
      Build: 8
      Sell:  3
      Hold:  5
Period: 2
  MultimediaBusiness:  
      Build: 11
      Sell:  16
      Hold:  0
  MultimediaHome:  
      Build: 9
      Sell:  7
      Hold:  5
  BasicHome:  
      Build: 0
      Sell:  4
      Hold:  1
Period: 3
  MultimediaBusiness:  
      Build: 17
      Sell:  17
      Hold:  0
  MultimediaHome:  
      Build: 4
      Sell:  9
      Hold:  0
  BasicHome:  
      Build: 4
      Sell:  5
      Hold:  0

```

表にしてまとめると以下になります。各期の組立台数、販売台数、在庫台数を振り分けました。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/223763/9c904c23-3e02-1d26-c809-28781669ebbc.png)





# 完成OPL

製品サンプルとほぼ同じ内容ですが、日本語コメントを入れたり順序を入れ替えたり、添え字を入れ替えたりして読みやすくしています。

```Java:totalprod_hkwd.mod
/*  ----------------------------------------------------
 *   OPL Model for Production planning Example
 *
 * This model is described in the documentation. 
 * See IDE and OPL > Language and Interfaces Examples.
　* https://www.ibm.com/docs/ja/icos/20.1.0?topic=library-manufacturing#usroplsamples.uss_opl_modlib.1016611__usroplsamples.uss_opl_modlib.1016871

このモデルは、コンピューター組み立て工場の生産計画を説明するものです。この工場は、MultimediaBusiness、MultimediaHome、および BasicHome という 3 つのコンピューター・モデルを組み立てることができます。
組み立ては、CPU やモデムなどのコンポーネントを適切なコンピューターに取り付けることで構成されます。このモデルでは、以下の 2 つのステップを示しています。

1. 組み立てる各タイプのコンピューターの総数と、必要なコンポーネントを取得する場所を決定します。(totalprod_hkwd.mod)
2. 組み立てるコンピューターの数を決定したら、各期間に組み立てる数、販売する数、および在庫として保持する数を決定します。(period_hkwd.mod)
 *   --------------------------------------------------- */
 
// コンピュータータイプとコンポーネントの定義
{string} ComputerTypes = ...;
{string} ComponentTypes = ...;

//int NumComputerTypes = card(ComputerTypes);

/*   ---------------------------------------------------
 *   Each computer type has a list of components, a 
 *   selling price and a maximum on the number allowed 
 *   to be held over in inventory to the next period.
 各コンピュータタイプには、コンポーネントのリスト、販売価格、および在庫可能最大数があります。
 *   --------------------------------------------------- */
tuple computersToBuild
{
  {string} components;
  int      price;
  int      maxInventory;
}

computersToBuild Computers[ComputerTypes] = ...;

// Number of Periods to use for this problem
// 期間数
int NbPeriods = ...;
range Periods = 1..NbPeriods;



/*   ---------------------------------------------------
 *   Each computer type has a maximum and minimum amount
 *   that can be sold in each period.  These values are
 *   used to calculate a range on the total to build of
 *   each computer type.  There is also a plant capacity
 *   for each period.
 各コンピューターの種類には、各期間に販売できる最大数量と最小数量があります。 これらの値は、各コンピューターの種類の合計構築範囲を計算するために使用されます。 各期間のプラント容量もあります。
 *   --------------------------------------------------- */
//期ごとの最大の組立可能数
int MaxBuildPerPeriod[1..NbPeriods] = ...;

//コンピュータの種類と期ごとの最小最大の需要量
int MinDemand[ComputerTypes][1..NbPeriods] = ...;
int MaxDemand[ComputerTypes][1..NbPeriods] = ...;



//全期分の最大組立数
int MaxBuild = 0;
//コンピュータの種類ごと全期分の最小最大の需要量
int TotalMaxDemand[ComputerTypes];
int TotalMinDemand[ComputerTypes];
//未使用float TotalBuild[ComputerTypes];
//for文で初期化
execute INITIALIZE {
  for(var p in Periods) {
    MaxBuild = MaxBuild + MaxBuildPerPeriod[p];
    for(var c in ComputerTypes) {
      TotalMaxDemand[c] = MaxDemand[c][p] + TotalMaxDemand[c];
      TotalMinDemand[c] = MinDemand[c][p] + TotalMinDemand[c];
    }
  }
}

//部品サプライヤー
{string} Suppliers = ...;
//サプライヤーのコスト変化ポイントの定義
int NbCostChanges = ...;         // 区分線形のブレイク数　Number of "breaks" in PLF
int SupplyCostCostChanges[1..NbCostChanges] = ...;
//未使用range CostChangespp = 1..NbCostChanges+1;

//コンポーネントのタイプ毎、サプライヤーごとのコストの傾斜
tuple componentInfo
{
  string supplier;
  int    costSlope[1..NbCostChanges+1];   //CostChangesppと同じであるべき。should be CostChangespp
}
{componentInfo} Components[ComponentTypes] = ...;

tuple supplierMatch {
  string         component;
  componentInfo  componentInformation;
}
{supplierMatch} ComponentSupplierMatches = { <i,j> | i in ComponentTypes, j in Components[i] };


////決定変数
/*  ----------------------------------------------------
 *   Variables:
 *   build --   How many of each computer type to build.
 *   NecessaryComponents -- Based on build
 *   SuppliedComponents -- How many of each component
 *         type are supplied by which supplier
 *   inHouse -- How many of each component type are
 *         manufactured in House.
 *   grossProfit -- linear function of build
 *   cost -- function of SuppliedComponents and inHouse
 *   --------------------------------------------------- */

//組立数 -- 各タイプのコンピューターを何台組み立てるか。
//dvar float+ Build[ComputerTypes];
dvar int+ Build[ComputerTypes];
//必要な各コンポーネント数 -- Buildに基づく 
//dvar float+ NecessaryComponents[ComponentTypes];
dvar int+ NecessaryComponents[ComponentTypes];
//コンポーネント供給数 -- 各コンポーネント タイプがどのサプライヤーによっていくつ供給されるか。
//dvar float+ SuppliedComponents[ComponentSupplierMatches];
dvar int+ SuppliedComponents[ComponentSupplierMatches];
/*
//内製コンポーネント数 -- 各コンポーネント タイプのいくつを社内で製造するか。
dvar float+ InHouse[ComponentTypes];
*/

////決定式
//コスト -- SuppliedComponents と inHouse 
dexpr float Cost =       
  sum(m in ComponentSupplierMatches)
    piecewise(i in 1..NbCostChanges) {
      m.componentInformation.costSlope[i] -> SupplyCostCostChanges[i];   
      m.componentInformation.costSlope[NbCostChanges+1]
    } SuppliedComponents[m];
    
//売上 -- ビルドの線形関数
dexpr float Sales=
  sum(p in ComputerTypes) Computers[p].price * Build[p];

//粗利益
dexpr float GrossProfit = Sales - Cost;

////目的関数 粗利益の最大化
//minimize Cost-GrossProfit;

maximize GrossProfit;

////制約
/*  ----------------------------------------------------
 *   constraints
 *   --------------------------------------------------- */
subject to
{
  //制約1: 各コンピューターは最大組み立て数以下
  ctPlantCapacity:
    sum(p in ComputerTypes) Build[p] <= MaxBuild;
      
  //制約2:各コンピューターは最大需要数以下
  //制約3:各コンピューターは最小需要数以上
  forall(p in ComputerTypes){
    ctComputerTypeMaxDemand: Build[p] <= TotalMaxDemand[p];
    ctComputerTypeMinDemand: Build[p] >= TotalMinDemand[p];      
  }  


  // Get the necessary components
  //制約4:組み立てるコンピュータータイプの台数と必要なコンポーネント数が一致する
/*
  forall(c in ComponentTypes)
    ctDetermineAmtNecessary:      
      sum(p in ComputerTypes: c in Computers[p].components) 
        Build[p] == NecessaryComponents[c];
*/
  forall(c in ComponentTypes)
    ctDetermineAmtNecessary:      
      sum(p in ComputerTypes: c in Computers[p].components) 
        Build[p] == NecessaryComponents[c];
         
 //制約5:供給されるコンポーネント数と必要なコンポーネント数が一致する
  forall(c in ComponentTypes)
    ctDetermineAmtSupplied:      
      sum(m in ComponentSupplierMatches: c == m.component) 
         SuppliedComponents[m] == NecessaryComponents[c];
 
/*
//自社供給されるコンポーネント数とコンポーネント供給数が一致する
  forall(m in ComponentSupplierMatches: m.componentInformation.supplier == "InHouse")
    ctMadeInHouse:      
      InHouse[m.component] == SuppliedComponents[m];
*/         
}

main {

/*   ---------------------------------------------------
 *   Solve the first model.  This piecewise-linear 
 *   program determines the total number of each computer 
 *   type to build, as well as from where to acquire the 
 *   necessary components.  The objective is to 
 *   maximize profit (really, minimize -profit).
 *   The quantity.build var values are stored to use in
 *   the second step.
 最初のモデルを解きます。 この区分線形プログラムは、構築する各コンピューター タイプの合計数と、必要なコンポーネントをどこから入手するかを決定します。 
 目的は、利益を最大化することです (実際には、利益を最小化することです)。
　quantity.totalbuild 変数値は、2 番目のステップで使用するために保存されます。
 *   --------------------------------------------------- */

   //数量決定モデルの設定
   var quantity = thisOplModel;
   quantity.generate();
   if (cplex.solve()) {

     var c;   
     var m; 
     var totalbuild = new Array();

     //粗利益出力
     writeln("Gross profit: ", cplex.getObjValue());
     // Output & get total to build of each type
     // 各コンピュータータイプ毎の組立台数の合計を算出して出力
     for(c in quantity.ComputerTypes) {
       totalbuild[c] = quantity.Build[c].solutionValue;
       writeln("   ", c, ":  ", totalbuild[c]);
     }
     //内製で生産するコンポーネント
     writeln("Components to build in-house:");
     for(m in quantity.SuppliedComponents){
       if(quantity.SuppliedComponents[m] > 0 
       && m.componentInformation.supplier == "InHouse") {
         writeln("  ", m.component, ": ", quantity.SuppliedComponents[m].solutionValue);
       }
    }      

    /*
     for(c in quantity.ComponentTypes){
       if(quantity.InHouse[c] > 0) {
         writeln("  ", c, ": ", quantity.InHouse[c].solutionValue);
       }
    }      
    */

     /*   ---------------------------------------------------
      *   Solve the second model.  This linear program
      *   determines the number of each computer type to
      *   build, sell and hold in each period.  The 
      *   objective is to find a feasible solution.
      2 番目のモデルを解きます。 この線形プログラムは、各期間に製造、販売、保持する各タイプのコンピューターの数を決定します。 
      目的は、実現可能な解決策を見つけることです。
      *   --------------------------------------------------- */

     // 期間最適化のモデル読込
     var source = new IloOplModelSource("period_hkwd.mod");
     var def = new IloOplModelDefinition(source);
     var newCplex = new IloCplex();
     var dates = new IloOplModel(def, newCplex);
     var data = new IloOplDataSource("period_hkwd.dat");
     dates.addDataSource(data);
     data = dates.dataElements;
     //モデル１の最適化結果でデータを書き換える。
     for(c in dates.ComputerTypes) {
       data.TotalBuild[c] = totalbuild[c];
     }
     dates.generate();
     if (newCplex.solve()) {
       // 出力　Output
       for(var p in dates.Periods) {
         writeln("Period: ", p)
         for(c in dates.ComputerTypes) {
           writeln("  ", c, ":  ")
           writeln("      Build: ", dates.Build[c][p].solutionValue)
           writeln("      Sell:  ", dates.Sell[c][p].solutionValue)
           writeln("      Hold:  ", dates.InStockAtEndOfPeriod[c][p].solutionValue)
         }
       }
     }
     dates.end();
   } else {
     writeln("Could not determine the total number of each computer type to build");
   }
}

```

データは以下です。

```Java:production_hkwd.dat
ComputerTypes = { MultimediaBusiness, MultimediaHome, BasicHome };
ComponentTypes = { CPU150, CPU250, Modem288, Modem56K, SoundCard, VideoCard, CDRW };

NbPeriods = 3;
//MaxBuildPerPeriod = [50 10 20];
MaxBuildPerPeriod = [35 20 25];

MinDemand = [[5,6,7],
         [4,5,6],
         [3,4,5]];

MaxDemand = [[15,16,17],
         [14,15,16],
           [13,14,15]];

Computers = [<{CPU250, Modem56K,SoundCard,VideoCard,CDRW},8000,5>,
                     <{CPU250,Modem56K,SoundCard,VideoCard},4000,5>,
                     <{CPU150,Modem288,VideoCard},2000,5>];

Suppliers = { CheapComponents ModemsRUs DiscoutCPUs CDsForASong InHouse };

NbCostChanges = 1;
SupplyCostCostChanges = [20];
          
Components = [ {<CheapComponents, [500,450]>, <DiscoutCPUs, [475,450]>},
          {<CheapComponents, [325,300]>, <DiscoutCPUs, [400,250]>},
          {<CheapComponents, [125,100]>, <ModemsRUs, [150,75]>},
          {<CheapComponents, [90,90]>, <ModemsRUs, [95,85]>},
          {<CheapComponents,  [50,50]>, <InHouse, [55,45]>},
          {<CheapComponents,  [75,75]>, <InHouse, [80,50]>},
          {<CheapComponents,  [105,90]>, <CDsForASong, [100,100]>}];

};
```

# 参考
#### 製造：生産計画- IBM Documentation 

https://www.ibm.com/docs/ja/icos/20.1.0?topic=library-manufacturing#usroplsamples.uss_opl_modlib.1016611__usroplsamples.uss_opl_modlib.1016871


#### CPLEXサンプルを読み解く記事一覧

https://qiita.com/kawada2017/items/bdb85f66e5f1de6039d6