À propos du problème du voyageur de commerce

Récemment, j'ai eu la chance d'entrer en contact avec «l'optimisation mathématique» grâce à mon travail. C'est une bonne idée, alors j'aimerais publier les résultats d'une petite touche privée.

Ce que j'ai fait était un exemple d'optimisation mathématique, le «problème du vendeur circulaire» Comme cela peut être mis en œuvre avec la pulpe d'environnement libre, j'ai essayé de l'implémenter tout en préparant l'environnement.

Les problèmes que je veux essayer sont les suivants.

◆ 2.12 Problème de vendeur de patrouille https://www.msi.co.jp/nuopt/docs/v19/examples/html/02-12-00.html

J'ai dessiné une image simple ci-dessous. A, B, C, D sont les points sur la carte.

map.png

Quel itinéraire a la plus petite distance totale en quittant A et en revenant d'un seul coup? Il semble que le problème d'optimisation s'appelle le «problème du vendeur de patrouille». Cela s'appelle la difficulté NP, et cela semble assez difficile.

Cependant, l'approche avec le modeleur appelé pulp qui est touché dans l'environnement python est gratuite !!! C'est facile à faire, je vais donc me référer à l'article suivant.

◆ Optimisation mathématique à partir du problème du voyageur de commerce https://qiita.com/panchovie/items/6509fb54e3d53f4766aa

La grande chose est que nous avons tout, de l'installation de pulp à l'échantillon de code. Mettons immédiatement le premier paramètre dans l'exemple de code. Qu'en est-il du code suivant qui semble compiler pour le moment?

import pulp as pp

#Définition du modèle d'optimisation
mip_model = pp.LpProblem("tsp_mip", pp.LpMinimize)

N = 4
BigM = 10000

c_ = [
      [0,6,5,5], # pointA > [pointA,pointB,pointC,pointD]
      [6,0,7,4], # pointB > [pointA,pointB,pointC,pointD]
      [5,7,0,3], # pointC > [pointA,pointB,pointC,pointD]
      [5,4,3,0], # pointD > [pointA,pointB,pointC,pointD]
      ]

#Définition variable
x = [[pp.LpVariable("x(%s,%s)"%(i, j), cat="Binary") for i in range(N)] for j in range(N)]
u = [pp.LpVariable("u(%s)"%(i), cat="Continuous", lowBound=1.0, upBound=(N)) for i in range(N)]

#Définition et enregistrement de l'indice d'évaluation (équation (1))
objective = pp.lpSum(c_[i][j] * x[i][j] for i in range(N) for j in range(N) if i != j)
mip_model += objective

#Expression conditionnelle(2)Enregistrement de
for i in range(N):
    mip_model += pp.lpSum(x[i][j] for j in range(N) if i != j) == 1

#Expression conditionnelle(3)Enregistrement de
for i in range(N):
    mip_model += pp.lpSum(x[j][i] for j in range(N) if i != j) == 1

for i in range(N):
    mip_model += x[i][i] == 0
    
#Expression conditionnelle(4) (Contrainte MTZ)
for i in range(N):
    for j in range(N):
        if i != j:
            mip_model += u[i] + 1.0 - BigM * (1.0 - x[i][j]) <= u[j]


#Effectuer l'optimisation
status = mip_model.solve()

#Comprendre les résultats
print("Status: {}".format(pp.LpStatus[status]))
print("Optimal Value [a.u.]: {}".format(objective.value()))

for i in range(N):
    for j in range(N):
        if i != j:
            print("x[%d][%d]:%f" % (i,j,x[i][j].value()))

for i in range(len(u)):
    print("u[%d] %f" % (i,u[i].value()))

Si vous l'exécutez dans cet état ...

Status: Infeasible
Optimal Value [a.u.]: 18.0
x[0][1]:0.000000
x[0][2]:1.000000
x[0][3]:0.000000
x[1][0]:0.999800
x[1][2]:0.000000
x[1][3]:0.000200
x[2][0]:0.000200
x[2][1]:0.000000
x[2][3]:0.999800
x[3][0]:0.000000
x[3][1]:1.000000
x[3][2]:0.000000

Il semble qu'il n'ait pas été bien calculé. Il semble que les conditions de contrainte ne soient pas bien satisfaites, donc après divers examens, il semble que la formule de contrainte MTZ soit étrange en l'état. Si vous écrivez les quatre formules pour x à u dans l'état final, ce sera comme suit.

u[0] < u[2]\\
u[1] < u[0]\\
u[2] < u[3]\\
u[3] < u[1]

Quand il s'agit de satisfaire tout cela,

u[0] < u[2] < u[3] < u[1] < u[0] < ...

Il augmentera à l'infini, mais comme la plage effective de u est comprise entre 1 et 4, il semble qu'il y ait une contradiction. Par conséquent, en supprimant la restriction sur la variable qui retourne au point A à la fin, reconstruisons la formule comme suit.

u[0] < u[2]\\
u[2] < u[3]\\
u[3] < u[1]

Cela modifiera l'implémentation et changera la formule de condition MTZ en:

#Expression conditionnelle(4) (Contrainte MTZ)
for i in range(N):
    for j in range(1,N):
        if i != j:
            mip_model += u[i] + 1.0 - BigM * (1.0 - x[i][j]) <= u[j]

Le code source est résumé ci-dessous.

sample_route.py


import pulp as pp

#Définition du modèle d'optimisation
mip_model = pp.LpProblem("tsp_mip", pp.LpMinimize)

N = 4
BigM = 10000

c_ = [
      [0,6,5,5], # pointA > [pointA,pointB,pointC,pointD]
      [6,0,7,4], # pointB > [pointA,pointB,pointC,pointD]
      [5,7,0,3], # pointC > [pointA,pointB,pointC,pointD]
      [5,4,3,0], # pointD > [pointA,pointB,pointC,pointD]
      ]

#Définition variable
x = [[pp.LpVariable("x(%s,%s)"%(i, j), cat="Binary") for i in range(N)] for j in range(N)]
u = [pp.LpVariable("u(%s)"%(i), cat="Continuous", lowBound=1.0, upBound=(N)) for i in range(N)]

#Définition et enregistrement de l'indice d'évaluation (équation (1))
objective = pp.lpSum(c_[i][j] * x[i][j] for i in range(N) for j in range(N) if i != j)
mip_model += objective

#Expression conditionnelle(2)Enregistrement de
for i in range(N):
    mip_model += pp.lpSum(x[i][j] for j in range(N) if i != j) == 1

#Expression conditionnelle(3)Enregistrement de
for i in range(N):
    mip_model += pp.lpSum(x[j][i] for j in range(N) if i != j) == 1

for i in range(N):
    mip_model += x[i][i] == 0
    
#Expression conditionnelle(4) (Contrainte MTZ)
for i in range(N):
    for j in range(1,N):
        if i != j:
            mip_model += u[i] + 1.0 - BigM * (1.0 - x[i][j]) <= u[j]


#Effectuer l'optimisation
status = mip_model.solve()

#Comprendre les résultats
print("Status: {}".format(pp.LpStatus[status]))
print("Optimal Value [a.u.]: {}".format(objective.value()))

for i in range(N):
    for j in range(N):
        if i != j:
            print("x[%d][%d]:%f" % (i,j,x[i][j].value()))

for i in range(len(u)):
    print("u[%d] %f" % (i,u[i].value()))

Le résultat est le suivant.

Status: Optimal
Optimal Value [a.u.]: 18.0
x[0][1]:0.000000
x[0][2]:1.000000
x[0][3]:0.000000
x[1][0]:1.000000
x[1][2]:0.000000
x[1][3]:0.000000
x[2][0]:0.000000
x[2][1]:0.000000
x[2][3]:1.000000
x[3][0]:0.000000
x[3][1]:1.000000
x[3][2]:0.000000
u[0] 1.000000
u[1] 4.000000
u[2] 2.000000
u[3] 3.000000

Il a convergé en toute sécurité. Les résultats ne sont pas graphiques,

Point A → Point C → Point D → Point B → Point A

Semble être le plus court et la distance est de 18. Où le problème a été répertorié

Point A → Point B → Point D → Point C → Point A

Donc, après tout, la distance est de 18 et elle semble être la plus courte, donc il y a quelques réponses.

L'optimisation mathématique n'est très récente que dans des domaines que nous n'avons pas encore abordés. Il existe un environnement où vous pouvez vous sentir libre de l'essayer, alors je vais profiter de cette occasion pour étudier diverses choses.

Recommended Posts

À propos du problème du voyageur de commerce
À propos du problème du vendeur de patrouille commandé
Python: j'ai essayé le problème du voyageur de commerce
Résolvez le problème du voyageur de commerce avec OR-Tools
J'ai essayé de mettre en œuvre le problème du voyageur de commerce
Pensez au problème de changement minimum
Problème de méthode de gravure et de voyageur de commerce
À propos du test
À propos de la file d'attente
Essayez de résoudre le problème du voyageur de commerce avec un algorithme génétique (théorie)
Essayez de résoudre le problème du voyageur de commerce avec un algorithme génétique (résultat de l'exécution)
À propos de la fonction Déplier
À propos de la commande de service
Résolution du problème du voyageur de commerce avec l'algorithme génétique (GA) et sa bibliothèque (vcopt)
Résolvez le problème du voyageur de commerce asymétrique Python avec la méthode de branche et de liaison
À propos de la matrice de confusion
À propos du modèle de visiteur
Examiner le double problème
Résolvez le problème de Monty Hall
À propos du module Python venv
À propos de la fonction enumerate (python)
Une histoire sur la façon de traiter le problème CORS
À propos de la compréhension du lecteur en 3 points [...]
À propos des composants de Luigi
À propos des fonctionnalités de Python
Optimisation de combinaison - Problème typique - Problème de vendeur circulaire
J'ai essayé "Implémentation d'un algorithme génétique (GA) en python pour résoudre le problème du voyageur de commerce (TSP)"
Une histoire dans laquelle l'algorithme est arrivé à une conclusion ridicule en essayant de résoudre correctement le problème du voyageur de commerce
[Python] Qu'est-ce que @? (À propos des décorateurs)
À propos de la valeur de retour de pthread_mutex_init ()
À propos de la valeur de retour de l'histogramme.
À propos du type de base de Go
À propos de la limite supérieure de threads-max
À propos de l'option moyenne de sklearn.metrics.f1_score
À propos du comportement de yield_per de SqlAlchemy
Illustration des résultats du problème du sac à dos
À propos de la taille des points dans matplotlib
À propos de la liste de base des bases de Python
Pensez grossièrement à la fonction de perte
Résolvez un simple problème de voyageur de commerce à l'aide d'une machine Boltzmann avec recuit simulé
Mon ami semble faire du python, alors pensez au problème ~ fizzbuzz ~