#2 Google OR-toolsを使って線形計画問題に制約条件を付ける

# モチベーション
「あたらしい数理最適化 python言語とGurobiで解く」
を教科書として、数理最適化の基本となるところを、肌感覚で理解したい。そのために、

- 基礎的な内容
- その例題演習

を書くことで、自身の理解と、新たな仕事の確保を目指します。私の大本のモチベーションは[ここ](https://qiita.com/amazyra/items/6b23d56148534ef1916e)を参照してください。


# この記事を読むとできること

- Google OR-Toolで制約条件付きの簡単な問題を解けるようにする

を目指します。

# ステップ1: 写経
例の如く、[Python用のStart Guide](https://developers.google.com/optimization/introduction/python?hl=ja)から、写経をしていきます。以下、引用。

```python:test.py
from __future__ import print_function
from ortools.linear_solver import pywraplp


def main():
    solver = pywraplp.Solver(
        'SolveSimpleSystem', pywraplp.Solver.GLOP_LINEAR_PROGRAMMING)
    # Create the variables x and y.
    x = solver.NumVar(-solver.infinity(), solver.infinity(), 'x')
    y = solver.NumVar(-solver.infinity(), solver.infinity(), 'y')
    # Constraint 1: x + 2y <= 14.
    constraint1 = solver.Constraint(-solver.infinity(), 14)
    constraint1.SetCoefficient(x, 1)
    constraint1.SetCoefficient(y, 2)

    # Constraint 2: 3x - y >= 0.
    constraint2 = solver.Constraint(0, solver.infinity())
    constraint2.SetCoefficient(x, 3)
    constraint2.SetCoefficient(y, -1)

    # Constraint 3: x - y <= 2.
    constraint3 = solver.Constraint(-solver.infinity(), 2)
    constraint3.SetCoefficient(x, 1)
    constraint3.SetCoefficient(y, -1)
    objective = solver.Objective()
    objective.SetCoefficient(x, 1)
    objective.SetCoefficient(y, 1)
    objective.SetMaximization()
    # Call the solver and display the results.
    solver.Solve()
    print('Solution:')
    print('x = ', x.solution_value())
    print('y = ', y.solution_value())

if __name__ == '__main__':
    main()
```
今回も、実行すればちゃんと最後まで通りました。ありがてぇ。

# ステップ2: 意味を理解する
さて、写経は済んだので意味を理解していきましょう。



## ソルバのインスタンス
```python

  solver = pywraplp.Solver(
        'SolveSimpleSystem', pywraplp.Solver.GLOP_LINEAR_PROGRAMMING)
  
```
これは明らかにソルバのインスタンスを立てていますね。前回同様、第一引数は、名前です。前回は調べなかった、`pywraplp.Solver.GLOP_LINEAR_PROGRAMMING`は何を意味するんでしょうか。GLOPについて調べたところ、[こんなページ](https://developers.google.com/optimization/lp/glop)を見つけました。どうやら、Googleが作った線形ソルバのようです。計算量も空間量も優れていると書かれていますが、何がどうすごいのかまでは書いていませんね...。とはいえ、このページの最後に`The Stigler diet`なる線形計画問題を見つけました。どうやら[このおっさん](https://ja.wikipedia.org/wiki/%E3%82%B8%E3%83%A7%E3%83%BC%E3%82%B8%E3%83%BB%E3%82%B9%E3%83%86%E3%82%A3%E3%82%B0%E3%83%A9%E3%83%BC)が提唱した、人が最低限の栄養を確保した際にかかるコストの最小値を問う問題のようです。

## 変数
```python
   # Create the variables x and y.
    x = solver.NumVar(-solver.infinity(), solver.infinity(), 'x')
    y = solver.NumVar(-solver.infinity(), solver.infinity(), 'y')
 
```
変数の設定。前回と違って、無限の定義域が設定されています。Cでよくある持ち方だと思いますが、`solver.infinty()`が無限を意味するようですね。
## 制約条件
```python
    # Constraint 1: x + 2y <= 14.
    constraint1 = solver.Constraint(-solver.infinity(), 14)
    constraint1.SetCoefficient(x, 1)
    constraint1.SetCoefficient(y, 2)

```
これが制約条件ですね（`Constraint`は制約を意味します）。制約条件を作る際に、その条件の値域を設定して、そこに係数付きで変数を放り込むようです。今回は不等式ですが、等式の場合はどうするんでしょうね。

```python
    # Constraint 1: x + 2y = 14.
    constraint1 = solver.Constraint(14, 14)
    constraint1.SetCoefficient(x, 1)
    constraint1.SetCoefficient(y, 2)
```

とすると、誤差が引っかかって制約条件を満たさない気がします。試してみましょう。まずはそのままの場合。

```python
Solution:
x =  5.999999999999998
y =  3.9999999999999996
>>>
```
明らかに誤差がついてきています。これは期待大。上式にしてやってみます。

```python
Solution:
x =  6.0
y =  4.0
>>>
```
やべぇ、ちゃんと解いた...。等式にすることでむしろ精度があがってます。制約がキツいせいで、逆に大域的最適解に素早く収束したんでしょうか。もうちょっときつくしてみましょうか。

```python
    # Constraint 3: x - y <= 2.
    constraint3 = solver.Constraint(-solver.infinity(), 2)
    constraint3.SetCoefficient(x, 1)
    constraint3.SetCoefficient(y, -1)
```
を、

```python
    # Constraint 3: x - y = 2.
    constraint3 = solver.Constraint(2, 2)
    constraint3.SetCoefficient(x, 1)
    constraint3.SetCoefficient(y, -1)
```
に変えてみましょう。完全に連立方程式を解かせていますが、あんまりきにしない。

```python
Solution:
x =  6.0
y =  4.0
>>>
```
ふむ。それでこそソルバだという気分にさせられます。よくできていて、面白い。しかし、なんで条件をきつくすると精度も上がるんですかね....。等式という強い制約条件を満たせるから（実行可能解があるから）さらに正確に解いた、が正解でしょうか。

なお、等式の取扱についてgoogle先生に聞いたところ、[Constraint Optimization](https://developers.google.com/optimization/cp/)なるものを教えてくれました。これは目的関数云々より、制約条件を満たす変数の組を探索してくれるもの。実行可能性を探してくれるともいえそうです。

## 目的関数

```python
    objective = solver.Objective()
    objective.SetCoefficient(x, 1)
    objective.SetCoefficient(y, 1)
    objective.SetMaximization()
```
ここは前と変わらず。大体わかってきました。

## 最適化の実行
```python
    solver.Solve()
    print('Solution:')
    print('x = ', x.solution_value())
    print('y = ', y.solution_value())
```
ここも疑問点はありませんね。pythonの、値渡しと参照渡しを同時に行う仕様がバリバリ使われていて、気を使うところでもありますね。値渡ししかしない世界に長いこと居すぎたので、下手に書くと偉いことになりそうです。

# まとめと考察
とりあえず、制約条件つきの線形計画問題は解くことができました。だいたいの構造も把握してきたので、そろそろ実践問題を一個解いてみたいですね。調べてる途中で行き着いた、以下の話はやって見る価値がありそうです。

- The Stigler diet
- Constraint Optimization
- その他、私が知らない最適化問題の形式


# 終わりに
自分の記事のクオリティを高めたいので、質問でも、次の要望、ただの感想、何でもいいのでコメントを残して頂けると幸いです。ちゃんと答えますので、よしなしによろしくお願いいたします。


# 関連記事
前回：
- [#1 Google OR-toolsを使って hello-worldする](https://qiita.com/amazyra/items/bcad1d23c83a6ad629b3)