diff --git a/exhaustive_search.py b/exhaustive_search.py
index 7c93d52..e770c2f 100644
--- a/exhaustive_search.py
+++ b/exhaustive_search.py
@@ -1,11 +1,12 @@
 import time
 from itertools import permutations
 from typing import Tuple
+from math import factorial
 
 import numpy as np
 import numpy.typing as npt
 
-from common import plot_plan, read_data
+from common import indexes_to_cities, plot_plan, read_data
 
 
 def exhaustive_search(distances: npt.NDArray) -> Tuple[float, npt.NDArray]:
@@ -34,14 +35,17 @@ def exhaustive_search(distances: npt.NDArray) -> Tuple[float, npt.NDArray]:
             ),
             permutations(range(size)),
         ),
+        key=lambda x: x[0] # Make sure that it finds the minimal distance
     )
 
 
 if __name__ == "__main__":
     cities, data = read_data("./european_cities.csv")
 
+    times = {}
+
     # A loop timing finding the optimal solution for different n
-    for n in range(6, 11):
+    for n in range(6,11):
         # Time exhaustive search
         t0 = time.time_ns()
         distance, perm = exhaustive_search(data[:n, :n])
@@ -49,32 +53,32 @@ if __name__ == "__main__":
 
         time_elapsed_ms = (t1 - t0) / 1_000_000.0
 
+        times[n] = time_elapsed_ms, distance
+
+        if n in (6,10):
+            city_seq = indexes_to_cities(perm, cities)
+            plot_plan(city_seq)
+            print(f"Sequence for {n} cities: {city_seq}")
+
+    print("")
+    for n, (time, distance) in times.items():
         print(f"Exhaustive search for the {n} first cities:")
-        print(f"distance             : {distance:>12.6f}km")
-        print(f"time to find solution: {time_elapsed_ms:>12.6f}ms\n")
+        print(f"{'distance':<25}: {distance:>12.6f}km")
+        print(f"{'time to find solution':<25}: {time:>12.6f}ms")
+        print(f"{f'time / {n}!':<25}: {time / factorial(n):>12.6f}\n")
 
 
 """Running example
 
-oblig1 on  main [?] via 🐍 v3.12.6 took 7s
+oblig1 on  main [!] via 🐍 v3.12.6 took 14s
 ❯ python exhaustive_search.py
 Exhaustive search for the 6 first cities:
-distance             :  5018.810000km
-time to find solution:     1.105330ms
-
-Exhaustive search for the 7 first cities:
-distance             :  5487.890000km
-time to find solution:    10.089604ms
-
-Exhaustive search for the 8 first cities:
-distance             :  6667.490000km
-time to find solution:    78.810508ms
-
-Exhaustive search for the 9 first cities:
-distance             :  6678.550000km
-time to find solution:   765.676230ms
+distance                 :  5018.810000km
+time to find solution    :     1.485208ms
+time / 6!                :     0.002063
 
 Exhaustive search for the 10 first cities:
-distance             :  7486.310000km
-time to find solution:  8281.795515ms
+distance                 :  7486.310000km
+time to find solution    : 10980.900480ms
+time / 10!               :     0.003026
 """
diff --git a/hill_climbing.py b/hill_climbing.py
index cc9cfda..6980af4 100644
--- a/hill_climbing.py
+++ b/hill_climbing.py
@@ -59,8 +59,29 @@ def hill_climbing(distances: npt.NDArray) -> Tuple[float, npt.NDArray]:
     return (current_distance, perm)
 
 
+def test_hill_climbing(data: npt.NDArray, cities: npt.NDArray, runs: int):
+    res = [hill_climbing(data) for _ in range(runs)]
+    res.sort(key=lambda n: n[0])
+
+    distances = list(map(lambda n: n[0], res))
+    best = res[0][0]
+    worst = res[-1][0]
+    avg = sum(distances) / runs
+    standard_deviation = np.sqrt(sum([(i - avg)**2 for i in distances]) / runs)
+
+    print(f"Hill climbing for {len(data)} cities.")
+    print(f"best distance     : {best:>12.6f}km")
+    print(f"worst distance    : {worst:>12.6f}km")
+    print(f"average distance  : {avg:>12.6f}km")
+    print(f"standard deviation: {standard_deviation:>12.6f}km\n")
+
+    plot_plan(indexes_to_cities(res[0][1], cities)) # Plot the best one
+
 if __name__ == "__main__":
+    np.random.seed(1987)
     cities, data = read_data("./european_cities.csv")
     distance, perm = hill_climbing(data[:10, :10])
 
-    plot_plan(indexes_to_cities(perm, cities))
+    # plot_plan(indexes_to_cities(perm, cities))
+    test_hill_climbing(data[:10,:10], cities, 20)
+    test_hill_climbing(data, cities, 20)