a3-algorithmique-avancee/Projet_algo.ipynb

Project Algorithmique avancée¶

Contexte du projet ¶

L'ADEME a lancé un appel pour promouvoir de nouvelles solutions de mobilité. CesiCDP, avec son équipe de 4 personnes, répond à cet appel en se concentrant sur la gestion des tournées de livraison. L'objectif est de minimiser la durée totale de la tournée en tenant compte des contraintes de trafic et autres. Un financement attractif est en jeu pour développer l'activité de CesiCDP. L'ajout de contraintes supplémentaires au problème rendra l'étude plus réaliste et plus difficile à résoudre, mais démontrera l'expertise de l'équipe.

Notre objectif est de développer un algorithme qui nous permettra de passer par tous les points de livraison dans un temps optimisé.

Contraintes choisies ¶

La contrainte que nous avons choisie est la suivante :

• Avoir plusieurs camions disponibles simultanément pour effectuer les livraisons.

Formulation du problème ¶

Considérons un graphe G=(V,E), où V est l'ensemble des villes (ou points de livraison) et E est l'ensemble des routes entre les villes. Il y a k camions disponibles pour effectuer les livraisons.

Le problème consiste à trouver un itinéraire pour chaque camion, de sorte que toutes les livraisons soient effectuées dans le temps le plus court possible, à la fois vers et depuis l'entrepôt.

Le problème que nous avons avec les contraintes mentionnées ci-dessus est le VRP (Problème de Routage des Véhicules).

Contraintes du problème ¶

Liste des contraintes du problème :

• Tous les clients doivent être servis
• Un client ne peut être servi que par un seul véhicule.
• En quittant un client, un véhicule ne peut se rendre que chez un seul autre client.

Nous allons donc affecter chaque client à une route desservie par un seul véhicule.

Démonstration de la complexité du problème de routage des véhicules (VRP) ¶

Introduction¶

Le problème du routage des véhicules (VRP) est une extension du problème du voyageur de commerce (TSP), et est connu pour être un problème NP-hard. Nous démontrerons la complexité du VRP en nous basant sur le problème de la chaîne hamiltonienne, qui est connu pour être NP-complet.

Définitions préliminaires¶

• Problème de la chaîne hamiltonienne** : le problème consiste à déterminer si un graphe non orienté donné possède une chaîne hamiltonienne, c'est-à-dire un chemin qui visite chaque sommet exactement une fois.

• Problème du voyageur de commerce (TSP)** : le problème consiste à trouver le circuit le plus court possible qui visite chaque ville d'un ensemble donné de villes et revient à la ville d'origine.

• Vehicle Routing Problem (VRP)** : le problème consiste à livrer une série de clients avec plusieurs véhicules tout en minimisant le coût total, tel que le temps de livraison ou la distance parcourue.

Preuve de la complexité du TSP.¶

Le TSP est une extension du problème de la chaîne hamiltonienne. En fait, un cas particulier du TSP est le problème de la chaîne hamiltonienne, dans lequel toutes les arêtes ont le même poids (ou coût). Dans ce cas, trouver le circuit le plus court dans le TSP est équivalent à trouver une chaîne hamiltonienne dans le graphe. Comme le problème de la chaîne hamiltonienne est NP-complet, le TSP doit être au moins aussi difficile, et le TSP est donc NP-complet.

Preuve de la complexité du VRP.¶

Le VRP est une extension du TSP. En fait, un cas particulier du VRP est le TSP où il n'y a qu'un seul véhicule disponible pour effectuer les livraisons. Dans ce cas, trouver la solution du VRP revient à trouver la solution du TSP.

Puisque le TSP est NP-complet, le VRP doit être au moins aussi difficile, donc le VRP est NP-difficile. En outre, le VRP introduit des contraintes supplémentaires, telles que des véhicules multiples et potentiellement des capacités de véhicules et des fenêtres temporelles, ce qui le rend encore plus complexe.

Conclusion¶

En conclusion, nous avons montré que le VRP est un problème NP-hard en le réduisant au problème NP-complet TSP, qui à son tour peut être réduit au problème NP-complet de la chaîne hamiltonienne. Cette démonstration met en évidence la difficulté de résoudre le VRP, en particulier pour les instances de grande taille. En pratique, des méthodes (méta)heuristiques et approximatives sont souvent utilisées pour résoudre le VRP, comme nous le verrons plus loin.

Modélisation mathématique ¶

Ensemble et paramètres¶

$V=\{0,1,2,...,n_v\}$ : l'ensemble des villes, où 0 est la base (ou le dépôt), $1,2,...,n_v$ sont les villes de livraison. $n_v+1$ sera le dépôt pour le retour.
$K=\{1,2,...,k\}$ : tous les camions.
$E$ représente les arcs entre deux clients $i,j \in V$
$G=(V,E)$ : le graphe avec V comme sommets (villes) et E comme arêtes
$d_{ij}$​ : distance (ou temps de déplacement) entre la ville i et la ville $j$ (coût du déplacement)
$M$ : une grande constante.

Variables de décision¶

$x_{ijk}$​ : binary variable worth 1 if truck $k$ moves from city $i$ to city $j$, and 0 otherwise.

Fonction objective¶

$$\min \sum_{k∈K} \sum_{⁡i∈V} \sum_{⁡j∈V} d_{ij} x_{ijk} $$

Contraintes VRP¶

• Chaque ville est visitée une fois et une seule : $$\sum_{⁡k \in K} \sum_{⁡j \in V} x_{ijk} = 1, \forall i \in V, i \ne 0$$

• Si un camion visite une ville, il doit la quitter : $$\sum_{i \in V} x_{ijk} = \sum_{j \in V} x_{ijk}, \forall k \in K, \forall i \in V, \forall j \in V $$
  

• Contrainte d'élimination des tournées partielles (pour garantir que chaque camion effectue une tournée complète) : $$\sum_{i \in S, j \notin S} x_{ijk} \geq 1, \forall k \in K, \forall \; sous-ensemble \; S \; de \; V, 0 \in S, S \ne V $$

Contraintes VRPTW¶

• Chaque ville est visitée une et une seule fois : $$\sum_{k \in K} \sum_{j \in V}x_{ijk} = 1, \forall i \in V, i \ne 0 $$

• Si un camion visite une ville, il doit la quitter : ∑𝑖∈𝑉𝑥𝑖𝑗𝑘=∑𝑗∈𝑉𝑥𝑖𝑗𝑘, ∀𝑘∈𝐾, ∀𝑖∈𝑉, ∀𝑗∈𝑉

• Contrainte d'élimination des sous-tours (pour s'assurer que chaque camion effectue un tour complet) : ∑𝑖∈𝑆,𝑗∉𝑆𝑥𝑖𝑗𝑘≥1, ∀𝑘∈𝐾, ∀𝑠𝑢𝑏𝑠𝑒𝑡𝑆𝑑𝑒𝑉, 0∈𝑆, 𝑆≠𝑉

• Chaque camion arrive dans la fenêtre temporelle de chaque ville : a_ik + t_ij ≤ a_jk + M(1 - x_ijk), ∀k∈K, ∀(i, j)∈A, ∀a_ik ∈ [e_i, l_i]

• Chaque camion quitte la ville après le temps de service : a_ik + s_i = d_ik, ∀k∈K, ∀i∈V, ∀a_ik ∈ [e_i, l_i]

Remarque : Ici, a_ik est l'heure d'arrivée du véhicule k au nœud i, t_ij est le temps de trajet du nœud i au nœud j, M est un grand nombre positif, x_ijk est une variable binaire qui vaut 1 si le véhicule k se déplace du nœud i au nœud j et 0 sinon, e_i est l'heure à laquelle un service peut commencer au plus tôt au nœud i, l_i est l'heure à laquelle un service peut commencer au plus tard au nœud i, et s_i est le temps de service au nœud i.

Ces contraintes tiennent compte non seulement de la nécessité de visiter chaque ville une fois, mais aussi de l'obligation d'arriver et de partir dans des fenêtres de temps spécifiques, ce qui rend le problème plus complexe que le VRP standard.

Resolution algorithm¶

Dans le cadre de notre projet, nous avons envisagé d'adopter trois algorithmes heuristiques distincts pour résoudre notre problème d'optimisation de trajet entre les villes. Notre objectif est de les comparer afin de déterminer celui qui nous fournira l'itinéraire le plus proche de l'optimal. Les trois algorithmes que nous allons mettre en œuvre sont les suivants :

• L'algorithme des colonies de fourmis
• L'algorithme du recuit simulé
• L'algorithme des essaims particulaires

Nous examinerons ensuite leur mise en œuvre et leur fonctionnement pour évaluer leur performance.

Algorithme par colonies de fourmis¶

L'algorithme par colonies de fourmis est une approche inspirée du comportement des fourmis dans la recherche de chemins optimaux. En utilisant des heuristiques et des phéromones, cet algorithme permet de résoudre des problèmes d'optimisation de manière efficace.

In [30]:

# Variables to edit:
nb_ville = 100
max_coords = 1000
nb_truck = 4
max_time = 4
nb_ants = 10
alpha = 1
beta = 6
evaporation = 0.5
intensification = 2

In [31]:

# Import necessary libraries
import matplotlib.pyplot as plt
import random, time
from tests.libs.clustering import split_tour_across_clusters
from tests.libs.aco import AntColony, total_distance

# Function to generate city coordinates
def generate_cities(nb, max_coords=1000):
    # Generate a list of nb cities, each with random coordinates within the range of max_coords
    return [random.sample(range(max_coords), 2) for _ in range(nb)]

# Calculate maximum time per cluster
max_time_per_cluster = max_time / nb_truck

# Track the time taken to generate cities
start_time_generate = time.time()
# Generate cities
cities = generate_cities(nb_ville, max_coords)
# Set the first city's coordinates at the center
cities[0] = [max_coords/2, max_coords/2]
stop_time_generate = time.time()

# Track the time taken to split cities across clusters
start_time_split = time.time()
# Split cities into clusters
clusters = split_tour_across_clusters(cities, nb_truck)
stop_time_split = time.time()

# Print out the time taken for each process
print("\n---- TIME ----")
print("Generate cities time: ", stop_time_generate - start_time_generate)
print("Split cities time: ", stop_time_split - start_time_split)

# Initialize a new figure for displaying paths
plt.figure()
plt.title("Ant Colony Optimization")
# List of colors for plotting
colors = [
    '#1f77b4',  # Bleu moyen
    '#ff7f0e',  # Orange
    '#2ca02c',  # Vert
    '#d62728',  # Rouge
    '#9467bd',  # Violet
    '#8c564b',  # Marron
    '#e377c2',  # Rose
    '#7f7f7f',  # Gris
    '#bcbd22',  # Vert olive
    '#17becf',  # Turquoise
    '#1b9e77',  # Vert Teal
    '#d95f02',  # Orange foncé
    '#7570b3',  # Violet moyen
    '#e7298a',  # Fuchsia
    '#66a61e',  # Vert pomme
    '#e6ab02',  # Jaune or
    '#a6761d',  # Bronze
    '#666666',  # Gris foncé
    '#f781bf',  # Rose clair
    '#999999',  # Gris moyen
]

# Initialize a list for storing the best routes
best_routes = []

print("\n---- CLUSTER ----")
# Iterate through each cluster
for i, cluster_indices in enumerate(clusters.values()):
    # Select a color for the cluster
    color = colors[i % len(colors)]
    
    # Retrieve city coordinates for the cluster
    cluster_cities = [cities[index] for index in cluster_indices]

    # Run the AntColony algorithm and store the best route
    ant_colony = AntColony(cluster_cities, n_ants=nb_ants, max_time=max_time_per_cluster, alpha=alpha, beta=beta, evaporation=evaporation, intensification=intensification)
    best_route = ant_colony.run()
    best_routes.append((best_route, color))

    # Print the total distance for the cluster
    print("Total distance for cluster {} ({} cities) : {}".format(i+1, len(cluster_cities), total_distance(best_route)))

# Calculate and print total distance for all clusters
full_total_distance = 0
for route, color in best_routes:
    full_total_distance += total_distance(route)
print("Total distance for all clusters: ", full_total_distance)

# Plot each route
for i, (route, color) in enumerate(best_routes):
    x = [city[0] for city in route]
    y = [city[1] for city in route]
    # Close the loop of the route
    x.append(x[0])
    y.append(y[0])
    # Plot the route
    plt.plot(x, y, color=color, marker='o', linestyle='-', label=f"Cluster {i}")

# Display the plot
plt.show()

---- TIME ----
Generate cities time:  0.0
Split cities time:  0.04548168182373047

---- CLUSTER ----
Total distance for cluster 1 (24 cities) : 2315.875564700583
Total distance for cluster 2 (27 cities) : 2161.8716667459307
Total distance for cluster 3 (28 cities) : 2309.9191397305735
Total distance for cluster 4 (28 cities) : 2261.7683336604277
Total distance for all clusters:  9049.434704837515

Algorithm du recuit simulé¶

L'algorithme du recuit simulé est une technique d'optimisation inspirée du processus de recuit utilisé dans la métallurgie. Il permet d'explorer de manière probabiliste l'espace des solutions en acceptant occasionnellement des mouvements qui peuvent initialement sembler moins favorables. Cette approche permet de trouver des solutions globales et d'éviter de rester bloqué dans des optima locaux.

In [32]:

# Variables to edit:
nb_ville = 100
max_coords = 1000
nb_truck = 4
temperature = 10000
cooling_rate = 0.9999
temperature_ok = 0.001

In [33]:

# Import necessary libraries
import matplotlib.pyplot as plt
import random, time
from tests.libs.clustering import split_tour_across_clusters
from tests.libs.simulated_annealing import SimulatedAnnealing, total_distance

# Function to generate city coordinates
def generate_cities(nb, max_coords=1000):
    # Generate a list of nb cities, each with random coordinates within the range of max_coords
    return [random.sample(range(max_coords), 2) for _ in range(nb)]

# Track the time taken to generate cities
start_time_generate = time.time()
# Generate cities
cities = generate_cities(nb_ville, max_coords)
# Set the first city's coordinates at the center
cities[0] = [max_coords/2, max_coords/2]
stop_time_generate = time.time()

# Track the time taken to split cities across clusters
start_time_split = time.time()
# Split cities into clusters
clusters = split_tour_across_clusters(cities, nb_truck)
stop_time_split = time.time()

# Print out the time taken for each process
print("\n---- TIME ----")
print("Generate cities time: ", stop_time_generate - start_time_generate)
print("Split cities time: ", stop_time_split - start_time_split)

# Initialize a new figure for displaying paths
plt.figure()
plt.title("Simulated Annealing")
# List of colors for plotting
colors = [
    '#1f77b4',  # Bleu moyen
    '#ff7f0e',  # Orange
    '#2ca02c',  # Vert
    '#d62728',  # Rouge
    '#9467bd',  # Violet
    '#8c564b',  # Marron
    '#e377c2',  # Rose
    '#7f7f7f',  # Gris
    '#bcbd22',  # Vert olive
    '#17becf',  # Turquoise
    '#1b9e77',  # Vert Teal
    '#d95f02',  # Orange foncé
    '#7570b3',  # Violet moyen
    '#e7298a',  # Fuchsia
    '#66a61e',  # Vert pomme
    '#e6ab02',  # Jaune or
    '#a6761d',  # Bronze
    '#666666',  # Gris foncé
    '#f781bf',  # Rose clair
    '#999999',  # Gris moyen
]

# Initialize a list for storing the best routes
best_routes = []

# Iterate through each cluster
for i, cluster_indices in enumerate(clusters.values()):
    # Select a color for the cluster
    color = colors[i % len(colors)]
    
    # Retrieve city coordinates for the cluster
    cluster_cities = [cities[index] for index in cluster_indices]

    # Run the SimulatedAnnealing algorithm and store the best route
    simulated_annealing = SimulatedAnnealing(cluster_cities, temperature=10000, cooling_rate=0.999, temperature_ok=0.01)
    best_route = simulated_annealing.run()
    best_routes.append((best_route, color))

    # Print the total distance for the cluster
    print("Total distance for cluster {} ({} cities) : {}".format(i+1, len(cluster_cities), total_distance(best_route)))

# Calculate and print total distance for all clusters
full_total_distance = 0
for route, color in best_routes:
    full_total_distance += total_distance(route)
print("Total distance for all clusters: ", full_total_distance)

# Plot each route
for i, (route, color) in enumerate(best_routes):
    x = [city[0] for city in route]
    y = [city[1] for city in route]
    # Close the loop of the route
    x.append(x[0])
    y.append(y[0])
    # Plot the route
    plt.plot(x, y, color=color, marker='o', linestyle='-', label=f"Cluster {i}")
# Add legend to the plot
plt.legend(loc="best")

# Display the plot
plt.show()

---- TIME ----
Generate cities time:  0.0009975433349609375
Split cities time:  0.04487872123718262
Total distance for cluster 1 (27 cities) : 2664.6125913825863
Total distance for cluster 2 (27 cities) : 2419.498935358293
Total distance for cluster 3 (27 cities) : 2792.878477974977
Total distance for cluster 4 (27 cities) : 2914.2725993190484
Total distance for all clusters:  10791.262604034904

Algorithme par essaim de particules¶

L'algorithme des essaims particulaires est une méthode d'optimisation basée sur le comportement collectif des essaims d'animaux, tels que les oiseaux ou les poissons. Il utilise des particules qui se déplacent dans l'espace de recherche et communiquent entre elles pour trouver des solutions optimales. Cet algorithme est efficace pour explorer l'espace des solutions et trouver des optima locaux et globaux.

In [34]:

# Variables to edit:
num_cities = 50
num_trucks = 5
num_particles = 50
max_iterations = 300
inertia_weight = 2
cognitive_weight = 1
social_weight = 1.5
incompleteness_probability = 0.8
start_time = 8
end_time = 18
infinite_distance_value = 1e6

In [49]:

# Importing necessary libraries
import random
import numpy as np
import networkx as nx
import matplotlib.pyplot as plt
import matplotlib
from IPython.display import display, Markdown
import matplotlib.cm as cm
from tests.libs.pso import Particle

# Function to generate a matrix of distances between cities
def generate_distance_matrix(num_cities, incompleteness_probability):
    # Initialize an empty matrix with zeros
    distances = np.zeros((num_cities, num_cities))
    # Iterate over each city
    for i in range(num_cities):
        # For each city, iterate over every other city to calculate distance
        for j in range(i+1, num_cities):
            # With certain probability, some distances will be missing (infinite)
            if random.random() < incompleteness_probability:
                distances[i, j] = float('inf')
                distances[j, i] = float('inf')
            else:
                # Assign random integer value to the distances between i and j
                distances[i, j] = random.randint(1, 10)
                distances[j, i] = distances[i, j]
    return distances

# Function to generate random time windows for each city
def generate_time_windows(num_cities, start_time, end_time):
    time_windows = []
    # For each city, generate a random start and end time
    for _ in range(num_cities):
        start = random.randint(start_time, end_time - 1)
        end = random.randint(start + 1, end_time)
        time_windows.append((start, end))
    return time_windows

# Function to generate a graph representation of the cities and the distances between them
def generate_city_graph(distances, time_windows):
    num_cities = distances.shape[0]
    # Create an empty graph
    G = nx.Graph()
    # Add cities as nodes with their corresponding time windows
    for i in range(num_cities):
        G.add_node(i+1, time_window=time_windows[i])
    # Add edges between cities with their corresponding distances
    for u in range(num_cities):
        for v in range(u+1, num_cities):
            if distances[u, v] != float('inf'):
                G.add_edge(u+1, v+1, distance=distances[u, v])
    return G

# Generate the distances and time windows
distances = generate_distance_matrix(num_cities, incompleteness_probability)
time_windows = generate_time_windows(num_cities, start_time, end_time)
# Generate the city graph
G = generate_city_graph(distances, time_windows)

# Initialize the particles and the global best position and cost
particles = []
global_best_position = np.random.permutation(range(1, num_cities+1))
global_best_cost = float('inf')
for _ in range(num_particles):
    particle = Particle(num_cities, num_trucks, distances, time_windows, infinite_distance_value)
    particles.append(particle)
    particle_cost = particle.evaluate_cost()
    if particle_cost < global_best_cost:
        global_best_cost = particle_cost
        global_best_position = particle.position.copy()

# Main loop for the particle swarm optimization algorithm
iteration = 0
best_costs = []
while iteration < max_iterations:
    for particle in particles:
        particle.update_velocity(global_best_position, inertia_weight, cognitive_weight, social_weight)
        particle.update_position()
        particle_cost = particle.evaluate_cost()
        if particle_cost < global_best_cost:
            global_best_cost = particle_cost
            global_best_position = particle.position.copy()
    best_costs.append(global_best_cost)
    iteration += 1

# Display the best solution found
best_position_str = ", ".join(str(node) for node in global_best_position)
best_cost_str = str(global_best_cost)
markdown_text = f"### Best solution found:\n\n- **City Positions**: {best_position_str}\n- **Total Cost**: {best_cost_str}"
display(Markdown(markdown_text))

# ---------------- Plot the best cost for each iteration ----------------------------------------------------------------------

# Compute the layout positions for the nodes using the spring layout algorithm
pos = nx.spring_layout(G)

# Create labels for each node in the graph
labels = {node: str(node) for node in G.nodes()}

# Draw the nodes of the graph using the computed positions
nx.draw_networkx_nodes(G, pos)

# Draw the edges of the graph using the computed positions
nx.draw_networkx_edges(G, pos)

# Draw the labels for the nodes using the computed positions and labels
nx.draw_networkx_labels(G, pos, labels)

# Set the title of the plot
plt.title("Graph of cities with time windows")

# Turn off the axis display
plt.axis("off")

# Display the plot
plt.show()

# ---------------- Plot Truck Paths ----------------------------------------------------------------------

# Initialize an empty list for truck paths
truck_paths = []

# Create a colormap for the trucks
cmap = matplotlib.colormaps['gist_ncar']
colors = [cmap(i/num_trucks) for i in range(num_trucks)]

# Split the global best position into separate paths for each truck
for truck in range(num_trucks):
    path = []
    for i in range(len(global_best_position)):
        if i % num_trucks == truck:
            city = global_best_position[i]
            path.append(int(city))
    path.append(path[0])  # Making the path circular by adding the first city to the end of the path
    truck_paths.append(path)

# Plot the truck paths on a graph
plt.figure(figsize=(10, 6))
plt.title("Truck Paths")
plt.axis("off")
for i in range(num_trucks):
    path = truck_paths[i]
    path_str = " -> ".join(str(city) for city in path)
    text = f"Truck {i+1}: {path_str}"
    display(Markdown(text))
    edges = [(path[j], path[j+1]) for j in range(len(path)-1)]
    nx.draw_networkx_edges(G, pos, edgelist=edges, edge_color=colors[i], label=f"Truck {i+1}")
    nx.draw_networkx_labels(G, pos, labels)
plt.legend()
plt.show()

# ---------------- Plot path cost for Each Truck ----------------------------------------------------------------------

# Initialize an empty list for the total cost of each path
path_costs = []

# Calculate the total cost of each path
for path in truck_paths:
    path_cost = 0
    for i in range(len(path) - 1):
        u, v = path[i] - 1, path[i + 1] - 1  # Adjust the indices
        # If the distance is infinite, use a replacement value.
        distance = distances[u, v] if distances[u, v] != float('inf') else infinite_distance_value
        path_cost += distance
    path_costs.append(path_cost)

# Plot a histogram of the path costs
plt.figure(figsize=(10, 6))
plt.bar(range(1, num_trucks + 1), path_costs, color=cm.rainbow(np.linspace(0, 1, num_trucks)))
plt.xlabel("Truck")
plt.ylabel("Path Cost")
plt.title("Path Cost for Each Truck")
plt.show()

# ---------------- Plot Evolution of the Best Cost Over Iterations ----------------------------------------------------------------------

# Plot a graph of the evolution of the best cost
plt.figure(figsize=(10, 6))
plt.plot(range(1, max_iterations + 1), best_costs)
plt.title("Evolution of the Best Cost Over Iterations")
plt.xlabel("Iteration")
plt.ylabel("Best Cost")
plt.grid(True)
plt.show()

Best solution found:¶

• City Positions: 22, 18, 5, 15, 27, 4, 36, 33, 16, 45, 17, 24, 50, 20, 48, 7, 1, 8, 3, 35, 40, 6, 25, 11, 43, 28, 34, 30, 38, 29, 26, 23, 13, 39, 47, 12, 2, 31, 37, 9, 49, 44, 19, 10, 41, 42, 32, 21, 46, 14
• Total Cost: 86000793.0

Truck 1: 22 -> 4 -> 17 -> 7 -> 40 -> 28 -> 26 -> 12 -> 49 -> 42 -> 22

Truck 2: 18 -> 36 -> 24 -> 1 -> 6 -> 34 -> 23 -> 2 -> 44 -> 32 -> 18

Truck 3: 5 -> 33 -> 50 -> 8 -> 25 -> 30 -> 13 -> 31 -> 19 -> 21 -> 5

Truck 4: 15 -> 16 -> 20 -> 3 -> 11 -> 38 -> 39 -> 37 -> 10 -> 46 -> 15

Truck 5: 27 -> 45 -> 48 -> 35 -> 43 -> 29 -> 47 -> 9 -> 41 -> 14 -> 27

Experimental design¶

Test de la variation du coût optimal en faisant varier le paramètre "inertia weight"¶

In [44]:

# List of inertia values to test
inertia_values = np.linspace(0.1, 5, 10) 

# Initialize list to store optimal cost for each inertia value
optimal_costs = []

# Loop over each inertia value
for inertia_weight in inertia_values:
    # Initialize particles and the global best cost
    particles = []
    global_best_position = np.random.permutation(range(1, num_cities+1))
    global_best_cost = float('inf')
    
    # Initialize particles and find the best cost
    for _ in range(num_particles):
        particle = Particle(num_cities, num_trucks, distances, time_windows, infinite_distance_value)
        particles.append(particle)
        particle_cost = particle.evaluate_cost()
        if particle_cost < global_best_cost:
            global_best_cost = particle_cost
            global_best_position = particle.position.copy()
    
    # Iterate to find the optimal solution
    iteration = 0
    best_costs = []
    while iteration < max_iterations:
        for particle in particles:
            particle.update_velocity(global_best_position, inertia_weight, cognitive_weight, social_weight)
            particle.update_position()
            particle_cost = particle.evaluate_cost()
            if particle_cost < global_best_cost:
                global_best_cost = particle_cost
                global_best_position = particle.position.copy()
        best_costs.append(global_best_cost)
        iteration += 1

    # Store the optimal cost for this inertia value
    optimal_costs.append(global_best_cost)

# Plot the evolution of the optimal cost as a function of inertia
plt.figure(figsize=(10, 6))
plt.plot(inertia_values, optimal_costs)
plt.title("Evolution of Optimal Cost as a Function of Inertia")
plt.xlabel("Inertia")
plt.ylabel("Optimal Cost")
plt.show()

Test de la variation du coût optimal en faisant varier le paramètre "cognitive weight"¶

In [45]:

# List of cognitive weight values to test
cognitive_values = np.linspace(0.1, 5, 10)  

# Initialize list to store optimal cost for each cognitive weight value
optimal_costs = []

# Loop over each cognitive weight value
for cognitive_weight in cognitive_values:
    # Initialize particles and the global best cost
    particles = []
    global_best_position = np.random.permutation(range(1, num_cities+1))
    global_best_cost = float('inf')
    
    # Initialize particles and find the best cost
    for _ in range(num_particles):
        particle = Particle(num_cities, num_trucks, distances, time_windows, infinite_distance_value)
        particles.append(particle)
        particle_cost = particle.evaluate_cost()
        if particle_cost < global_best_cost:
            global_best_cost = particle_cost
            global_best_position = particle.position.copy()
    
    # Iterate to find the optimal solution
    iteration = 0
    best_costs = []
    while iteration < max_iterations:
        for particle in particles:
            particle.update_velocity(global_best_position, inertia_weight, cognitive_weight, social_weight)
            particle.update_position()
            particle_cost = particle.evaluate_cost()
            if particle_cost < global_best_cost:
                global_best_cost = particle_cost
                global_best_position = particle.position.copy()
        best_costs.append(global_best_cost)
        iteration += 1

    # Store the optimal cost for this cognitive weight value
    optimal_costs.append(global_best_cost)

# Plot the evolution of the optimal cost as a function of cognitive weight
plt.figure(figsize=(10, 6))
plt.plot(cognitive_values, optimal_costs)
plt.title("Evolution of Optimal Cost as a Function of Cognitive Weight")
plt.xlabel("Cognitive Weight")
plt.ylabel("Optimal Cost")
plt.show()

La méthode du coude¶

In [5]:

from sklearn.cluster import KMeans
import matplotlib.pyplot as plt
import numpy as np
import random, time, math

# Function to generate random coordinates for cities
def generate_cities(nb, max_coords=1000):
    return [random.sample(range(max_coords), 2) for _ in range(nb)]

nb_ville = 100  # number of cities
max_coords = 1000  # maximum value for coordinates

max_clusters = 20  # maximum number of clusters

start_time_generate = time.time()
cities = generate_cities(nb_ville, max_coords)  # generate cities
stop_time_generate = time.time()

# Function to apply the elbow method for determining optimal number of clusters
def elbow_method(cities, max_clusters):
    inertias = []  # list to store the inertia for each number of clusters
    for i in range(1, max_clusters+1):
        # Apply KMeans clustering
        kmeans = KMeans(n_clusters=i, random_state=0).fit(cities)
        # Append the inertia (sum of squared distances to closest centroid) to the list
        inertias.append(kmeans.inertia_)
    return inertias

cities = generate_cities(nb_ville, max_coords)  # generate cities 
inertias = elbow_method(cities, max_clusters)  # apply elbow method

plt.figure()
plt.plot(range(1, len(inertias)+1), inertias, marker='o')  # plot the inertia for each number of clusters
plt.title('Elbow Method')  # set title of the plot
plt.xlabel('Number of clusters')  # set x-axis label
plt.ylabel('Inertia')  # set y-axis label
plt.show()  # display the plot

c:\Python310\lib\site-packages\sklearn\cluster\_kmeans.py:870: FutureWarning: The default value of `n_init` will change from 10 to 'auto' in 1.4. Set the value of `n_init` explicitly to suppress the warning
  warnings.warn(
c:\Python310\lib\site-packages\sklearn\cluster\_kmeans.py:870: FutureWarning: The default value of `n_init` will change from 10 to 'auto' in 1.4. Set the value of `n_init` explicitly to suppress the warning
  warnings.warn(
c:\Python310\lib\site-packages\sklearn\cluster\_kmeans.py:870: FutureWarning: The default value of `n_init` will change from 10 to 'auto' in 1.4. Set the value of `n_init` explicitly to suppress the warning
  warnings.warn(
c:\Python310\lib\site-packages\sklearn\cluster\_kmeans.py:870: FutureWarning: The default value of `n_init` will change from 10 to 'auto' in 1.4. Set the value of `n_init` explicitly to suppress the warning
  warnings.warn(
c:\Python310\lib\site-packages\sklearn\cluster\_kmeans.py:870: FutureWarning: The default value of `n_init` will change from 10 to 'auto' in 1.4. Set the value of `n_init` explicitly to suppress the warning
  warnings.warn(
c:\Python310\lib\site-packages\sklearn\cluster\_kmeans.py:870: FutureWarning: The default value of `n_init` will change from 10 to 'auto' in 1.4. Set the value of `n_init` explicitly to suppress the warning
  warnings.warn(
c:\Python310\lib\site-packages\sklearn\cluster\_kmeans.py:870: FutureWarning: The default value of `n_init` will change from 10 to 'auto' in 1.4. Set the value of `n_init` explicitly to suppress the warning
  warnings.warn(
c:\Python310\lib\site-packages\sklearn\cluster\_kmeans.py:870: FutureWarning: The default value of `n_init` will change from 10 to 'auto' in 1.4. Set the value of `n_init` explicitly to suppress the warning
  warnings.warn(
c:\Python310\lib\site-packages\sklearn\cluster\_kmeans.py:870: FutureWarning: The default value of `n_init` will change from 10 to 'auto' in 1.4. Set the value of `n_init` explicitly to suppress the warning
  warnings.warn(
c:\Python310\lib\site-packages\sklearn\cluster\_kmeans.py:870: FutureWarning: The default value of `n_init` will change from 10 to 'auto' in 1.4. Set the value of `n_init` explicitly to suppress the warning
  warnings.warn(
c:\Python310\lib\site-packages\sklearn\cluster\_kmeans.py:870: FutureWarning: The default value of `n_init` will change from 10 to 'auto' in 1.4. Set the value of `n_init` explicitly to suppress the warning
  warnings.warn(
c:\Python310\lib\site-packages\sklearn\cluster\_kmeans.py:870: FutureWarning: The default value of `n_init` will change from 10 to 'auto' in 1.4. Set the value of `n_init` explicitly to suppress the warning
  warnings.warn(
c:\Python310\lib\site-packages\sklearn\cluster\_kmeans.py:870: FutureWarning: The default value of `n_init` will change from 10 to 'auto' in 1.4. Set the value of `n_init` explicitly to suppress the warning
  warnings.warn(
c:\Python310\lib\site-packages\sklearn\cluster\_kmeans.py:870: FutureWarning: The default value of `n_init` will change from 10 to 'auto' in 1.4. Set the value of `n_init` explicitly to suppress the warning
  warnings.warn(
c:\Python310\lib\site-packages\sklearn\cluster\_kmeans.py:870: FutureWarning: The default value of `n_init` will change from 10 to 'auto' in 1.4. Set the value of `n_init` explicitly to suppress the warning
  warnings.warn(
c:\Python310\lib\site-packages\sklearn\cluster\_kmeans.py:870: FutureWarning: The default value of `n_init` will change from 10 to 'auto' in 1.4. Set the value of `n_init` explicitly to suppress the warning
  warnings.warn(
c:\Python310\lib\site-packages\sklearn\cluster\_kmeans.py:870: FutureWarning: The default value of `n_init` will change from 10 to 'auto' in 1.4. Set the value of `n_init` explicitly to suppress the warning
  warnings.warn(
c:\Python310\lib\site-packages\sklearn\cluster\_kmeans.py:870: FutureWarning: The default value of `n_init` will change from 10 to 'auto' in 1.4. Set the value of `n_init` explicitly to suppress the warning
  warnings.warn(
c:\Python310\lib\site-packages\sklearn\cluster\_kmeans.py:870: FutureWarning: The default value of `n_init` will change from 10 to 'auto' in 1.4. Set the value of `n_init` explicitly to suppress the warning
  warnings.warn(
c:\Python310\lib\site-packages\sklearn\cluster\_kmeans.py:870: FutureWarning: The default value of `n_init` will change from 10 to 'auto' in 1.4. Set the value of `n_init` explicitly to suppress the warning
  warnings.warn(