容量制約付き配送計画問題（capacitated vehicle routing problem,
CVRP）に関して、以下のURLを参考にして定式化およびpythonでのシミュレーションを試みました。

参考文献：配送計画問題をpythonで最適化する  
URL:<https://www.letsopt.com/entry/2020/08/30/180258>  
作成日：2020/8/30  
作成者：cresselia2012

しかし、試みる過程で、以下のMTZ制約に関する制約式とpythonのプログラムがなぜ一致しているのか理解できませんでした。  
参考にさせていただいたURL内のpythonのプログラム、特にMTZ制約が正しい理由をご教授いただけますと大変ありがたいです。

理解できなかった理由としては2点あります。

  1. 参考文献に記載されている本来のMTZ制約の制約式ではu[j]とdemand[j]に入力される顧客jは同様のはずですが、参考文献に記載されているpythonのプログラムではu[j-1]とdemand[j]と記述されており、一見異なる顧客を入力しているように思えます。

[![\[参考文献で使用されている定式化\]](https://i.stack.imgur.com/nub6I.jpg)](https://i.stack.imgur.com/nub6I.jpg)

```

    ### MTZ制約
    for i in range(1,num_nodes):
        for j in range(1,num_nodes):
            if i != j and demand[i] + demand[j] <= capacity:
                problem += u[i-1] - u[j-1] + capacity * x[i][j] <= capacity - demand[j]
    
```

  2. 参考文献に記載されている変数u[i]の入力iは0を含まないと定義されています。  
しかし、MTZ制約に関するpythonのプログラムではi=1のとき、u[0]が生み出されてしまい、矛盾が起きてしまうのかなと疑問に思いました。

```

    ## 変数を定義
    u = [ pulp.LpVariable( 'u_{}'.format( i ), demand[i], capacity, cat="Integer" ) \
            for i in range(1,num_nodes) ]
    
    ### MTZ制約
    for i in range(1,num_nodes):
        for j in range(1,num_nodes):
            if i != j and demand[i] + demand[j] <= capacity:
                problem += u[i-1] - u[j-1] + capacity * x[i][j] <= capacity - demand[j]
    
```

ちなみに、以下のように本来のMTZ制約の制約式のままでpythonのプログラムを記述して、シミュレーションをしたところ、

```

    for i in range(1,num_nodes):
        for j in range(1,num_nodes):
            if i != j and demand[i] + demand[j] <= capacity:
                problem += u[i] - u[j] + capacity * x[i][j] <= capacity - demand[j]
    
```

結果としては以下のようなエラーが起き、プログラムが正しくないことがわかりました。  
u[i-1] - u[j-1]をu[i] - u[j]に変更したことでエラーが起こったのは明らかなのですが、具体的な原因はわかりませんでした。

```

    Traceback (most recent call last):
      File "C:\Users\Control-PC\Documents\easy_vrp5_compver.py", line 62, in <module>
        problem += u[i] - u[j] + capacity * x[i][j] <= capacity - demand[j]
    IndexError: list index out of range
    
```

最後に私がURLを参考にして記述したpythonのプログラムを記載させていただきます。  
以下のプログラムのMTZ制約は参考文献と同様のプログラムを使用しており、シミュレーションを行うとエラーは起きません。

```

    import math
    
    def makeCVRP():
        num_nodes = 7
        capacity = 30
        demand = [ 0, 9, 11, 13, 7,\
                   19, 17 ]
        coordinate = []
        coordinate.append( ( 190, 3 ) )
        coordinate.append( ( 98, 184 ) )
        coordinate.append( ( 5, 42 )  )
        coordinate.append( ( 117, 89 )  )
        coordinate.append( ( 61, 162 ) )
        coordinate.append( ( 9, 97 )  )
        coordinate.append( ( 80, 15 ) )
        def computeDistance( c1, c2 ):
            return math.sqrt( pow( c2[0] - c1[0], 2 ) + pow( c2[1] - c1[1], 2 ) )
        distance = [ [ round(computeDistance( c1, c2 )) for c1 in coordinate ] \
                        for c2 in coordinate ]
        return num_nodes, capacity, demand, distance, coordinate
    
    # make problem
    num_nodes, capacity, demand, distance, coordinate = makeCVRP()
    
    import pulp
    
    
    # 最適化問題を定義
    problem = pulp.LpProblem( "CVRP", pulp.LpMinimize )
    
    ## 変数を定義
    x = [ [ pulp.LpVariable( 'x_{}_{}'.format( i, j ), cat="Binary" ) \
            if i != j else None for j in range(num_nodes) ] \
            for i in range(num_nodes) ]
    u = [ pulp.LpVariable( 'u_{}'.format( i ), demand[i], capacity, cat="Integer" ) \
            for i in range(1,num_nodes) ]
    num_vehicle= [ pulp.LpVariable( 'num_vehicle', 1, num_nodes, cat="Integer" ) ]
    
    ## 最小化したい関数を定義
    problem += pulp.lpSum( distance[i][j] * x[i][j] for i in range(num_nodes) \
                            for j in range(num_nodes) if i != j )
    
    ## 制約条件を定義
    ### \sum x_ij = 1
    for j in range(1,num_nodes):
        problem += pulp.lpSum( x[i][j] for i in range(num_nodes) if i != j ) == 1
    
    for i in range(1,num_nodes):
        problem += pulp.lpSum( x[i][j] for j in range(num_nodes) if i != j ) == 1
    
    ### \sum x_0j = |K|
    problem += pulp.lpSum(x[0][j] for j in range(1,num_nodes)) == num_vehicle
    
    ### MTZ制約
    for i in range(1,num_nodes):
        for j in range(1,num_nodes):
            if i != j:
                problem += u[i-1] - u[j-1] + capacity * x[i][j] <= capacity - demand[j]
    
    # solve
    result = problem.solve(pulp.CPLEX_CMD())
    
    print("objective value = {}".format(pulp.value(problem.objective)))
    
    
    # pulp（CBC）の結果から辺を定義
    edges = [ ( i, j ) for i in range(num_nodes) for j in range(num_nodes)
                if i != j and pulp.value(x[i][j]) == 1 ]
                
    
    # edgesから経路を取得
    paths = []
    for i,j in edges:
        if i == 0:
            path = [ i, j ]
            while path[-1] != 0:
                for v, u in edges:
                    if v == path[-1]:
                        path.append(u)
                        break
            paths.append(path)
    
    import networkx as nx
    import matplotlib.pyplot as plt
    
    # 有向グラフの作成
    G = nx.DiGraph()
    G.add_edges_from( edges )
    
    color = [ "r", "g", "y", "m", "c" ]
    edge_color = []
    
    # 経路毎に色をセット
    for i,j in G.edges:
        for t,path in enumerate(paths):
            if i in path and j in path:
                edge_color.append( color[t] )
                break
    assert len(edges) == len(edge_color)
    
    # グラフの描画
    pos = { i : coordinate[i] for i in range(num_nodes) }
    
    fig = plt.figure()
    nx.draw_networkx( G, pos, edge_color=edge_color, alpha=0.5)
    
    # 画像を保存
    plt.axis("off")
    fig.savefig("test.png")
    
    print(plt.show())
    
```

初めての質問であり、質問の仕方が間違っていたら申し訳ありません。  
何卒よろしくお願い申し上げます。
CC BY-SA 4.0
2021-07-19T08:54:10.660
78292
2021-07-20T01:34:01.090
47353
question
python
容量制約付き配送計画問題（CVRP）の定式化とpythonプログラムについて
大変有難いことに知見のある方からご教授いただきました。  
理解できていなかった原因を単刀直入に申し上げますと、私のpythonの基礎的な知識不足でした。(そのため専門性は必要ありませんでした。)

以下にご教授いただいた内容を段階的に記載いたします。

step1:  
uとdのpython内の定義より、uの添字0の場所(uの配列内の0番目)にu[1]の値が入っており、またdの添字1の場所(dの配列内の1番目)にd[1]の値が入っている。  
step2:  
よって、u[1]の値とd[1]の値を対応させるためには、uの添字0の場所(uの配列内の0番目)からu[1]の値を、dの添字1の場所からd[1]の値を取り出す必要がある。  
step3:  
そのため、事前に定式化したMTZ制約を記述した際はu[i] - u[j] + capacity * x[i][j] <= capacity -
demand[j]であったが、  
今回の場合(uとdのpython内の定義的に)、事前に定式化したMTZ制約と同様の意味を持つMTZ制約をpythonで記述するには、uの添字内を-1し、u[i-1]
- u[j-1] + capacity * x[i][j] <= capacity - demand[j]と記述する必要がある。

（例えばj=1をpythonのプログラムに入れてみると、uの添字0の場所から変数u[1]を、dの添字1の場所から変数d[1]を取り出すことができ、u[j-1]とdemand[j]が対応していることがわかる。）

最後に、長い質問文をお読み下さりましたユーザー様、ご教授下さりました回答者様に、御礼申し上げます。
CC BY-SA 4.0
2021-07-20T01:34:01.090
78307
2021-07-20T01:34:01.090
47353
78292
answer