Préparez ce dont vous avez besoin.

python

import random, numpy as np, pandas as pd, networkx as nx, matplotlib.pyplot as plt
from itertools import chain, combinations
from pulp import *

def draw(g):
    """dessin"""
    nx.draw_networkx_labels(g, pos=pos)
    nx.draw_networkx_nodes(g, node_color='w', pos=pos)
    nx.draw_networkx_edges(g, pos=pos)
    plt.show()
def addvar(cnt=[0], *args, **kwargs):
    """Création de variables"""
    cnt[0] += 1
    return LpVariable('v%d'%cnt[0], lowBound=0, *args, **kwargs)

Créons un problème au hasard. La demande n'a été établie qu'entre certains nœuds.

python

n = 16 #Nombre de nœuds
g = nx.random_graphs.fast_gnp_random_graph(n, 0.26, 8) #Graphique
rn = g.nodes() #Liste des nœuds
pos = nx.spring_layout(g, pos={i:(i/4, i%4) for i in rn}) #Position du nœud
for i, j in g.edges():
    v = pos[i] - pos[j]
    g[i][j]['dist'] = np.sqrt(v.dot(v)) #distance
dems = random.sample(list(combinations(rn, 2)), 10) #demande
draw(g)

J'essaierai de trouver l'arbre de la superficie totale minimale.

python

h = nx.minimum_spanning_tree(g, 'dist')
draw(h)

J'essaierai de le résoudre.

python

rved = [(j, i) for i, j in g.edges()] #Verso
#Débit par demande
a = pd.DataFrame([(df, dt, ef, et, g.edge[ef][et]['dist'], addvar()) 
                  for df, dt in dems for ef, et in chain(g.edges(), rved)], 
                 columns=['DeFr', 'DeTo', 'EdFr', 'EdTo', 'Dist', 'Var'])
#S'il faut installer un côté
b = pd.DataFrame([(fr, to, g.edge[fr][to]['dist'], addvar(cat=LpBinary)) 
                  for fr, to in g.edges()], columns=['Fr', 'To', 'Dist', 'Var'])
res = [] #Solution(coût d'installation,Commodité,Graphique)
mxcst = 999 #Limite supérieure du coût d'installation
while True:
    m = LpProblem() #Modèle mathématique
    m += lpDot(a.Dist, a.Var) + lpDot(b.Dist, b.Var)*1e-6 #Fonction objective(Commodité)
    m += lpDot(b.Dist, b.Var) <= mxcst #Limite supérieure du coût d'installation
    for _, r in a.iterrows():
        i, j = r.EdFr, r.EdTo
        if i > j: i, j = j, i
        #Installer lors de l'écoulement du flux
        m += r.Var <= lpSum(b.query('Fr==%s & To==%s'%(i,j)).Var)
    for (df, dt), c in a.groupby(('DeFr', 'DeTo')):
        for nd in rn:
            z = 1 if nd == df else -1 if nd == dt else 0
            #Stockage de flux
            m += lpSum(c.query('EdFr == %s'%nd).Var) == \
                 lpSum(c.query('EdTo == %s'%nd).Var) + z
    m.solve() #Solution
    if m.status != 1: break
    a['Val'] = a.Var.apply(lambda v: value(v)) #résultat(Débit)
    b['Val'] = b.Var.apply(lambda v: value(v)) #résultat(S'il faut installer)
    cst = value(lpDot(b.Dist, b.Var)) #coût d'installation
    val = value(m.objective) #Commodité
    mxcst = cst - 1e-3 #Limite supérieure du coût d'installation suivant
    h = nx.Graph()
    h.add_edges_from([(r.Fr, r.To) for _, r in b[b.Val == 1].iterrows()])
    res.append((cst, val, h))

Voyons le résultat. Plus le coût d'installation est élevé, meilleure est la commodité.

python

plt.rcParams['font.family'] = 'IPAexGothic'
plt.plot([r[0] for r in res], [r[1] for r in res])
plt.xlabel('coût d'installation')
plt.ylabel('Commodité')

C'est le graphique du cas le plus pratique.

python

draw(res[0][2])

J'ai essayé de faire une animation GIF de la transition du coût d'installation.

c'est tout