Added 'simu-opti' practical work.

simu-opti/simu.py

@@ -0,0 +1,107 @@
+import networkx as nx
+import matplotlib.pyplot as plt
+import numpy as np
+import itertools
+import math
+
+MAX = 10
+print("Launching with ", MAX, " iterations...")
+import numpy as np
+
+def write_to_file(P_sequences, filename="P_sequences.txt"):
+    with open(filename, "w") as file:
+        for idx, P_sequence in enumerate(P_sequences):
+            file.write(f"P_sequence {idx + 1}:\n")
+            np.savetxt(file, P_sequence, fmt="%.6f")
+            file.write("\n")  # Add a blank line between matrices
+
+def converge(v,p):
+    i = 0
+    while (np.max(np.abs(v - np.mean(v))) > 1):
+        v = np.matmul(p, v)
+        i += 1
+    return i
+
+def create_vector(max):
+    np.random.seed(0)  # For reproducibility
+    vector = np.random.randint(0, 1000, size=(max, 1))  # Random integers
+    while np.any(np.abs(np.diff(vector.flatten())) < 100):  # Ensure the difference is at least 100
+        vector = np.random.randint(0, 1000, size=(max, 1))
+    return vector
+
+# Creating graph
+G = nx.Graph()
+for i in range(MAX-1):
+    G.add_edge(i+1, i+2)
+
+# Initialize the matrix for each node
+E_mat = {}
+for i in G.nodes:
+    # 0 either way for each value of the matrix
+    E = np.zeros((MAX, MAX))
+
+
+    for n in G.nodes:
+        for p in G.nodes:
+            # for each node -> n, then each node -> p
+            # n is not a neightbor [v_n does not belong to N(v_i)]
+            if n not in G.neighbors(i): 
+                if p == n:
+                    E[n - 1, p - 1] = 1
+            # n is a neightbor [v_n belongs to N(v_i)]
+            elif n in G.neighbors(i):
+                if p == n or p == i:
+                    E[n - 1, p - 1] = 0.5
+
+    # Store the matrix for this node
+    E_mat[i] = E
+
+# Emission sequence
+# For example sigma = sigma_1, sigma_2, ..., sigma_N we have :
+
+# Creating multiple possible iterations
+node_permutations = itertools.permutations(range(1, MAX + 1))
+
+i = 0
+P_sequences = []
+
+# Eighenvalues high -> worst scenario, Eighenvalues low -> best scenario (lowest time)
+lambda_2_values = []
+convergence_times = []
+
+for perm in node_permutations:
+    # calculate sequence
+    P_sequence = E_mat[perm[0]]
+    for j in range(1,MAX):
+        P_sequence = np.matmul(P_sequence, E_mat[perm[j]])
+    
+    # calculate eighenvalues and convergence time
+    eigenvalues, _ = np.linalg.eig(P_sequence)
+    # sort them then get best second value
+    sorted_eigenvalues = np.sort(np.abs(eigenvalues))[::-1]
+    lambda_2 = sorted_eigenvalues[1]
+    # storing them by calculating the convergence time (1/1-lambda_2)
+    lambda_2_values.append({'p_seq': P_sequence, 'val': lambda_2})
+    #convergence_times.append(converge(vector_0,P_sequence))
+
+    i = i+1
+
+print("Theorical: ", math.factorial(MAX))
+print("Final i: ",i)
+
+# Create the vector v0
+vector_0 = create_vector(MAX)
+
+lambda_2_values.sort(key=lambda x: x['val'])
+print("\n\nBest case: ", converge(vector_0,lambda_2_values[0]['p_seq']))
+print("Worst case: ", converge(vector_0,lambda_2_values[len(lambda_2_values)-1]['p_seq']))
+#write_to_file(P_sequences, filename="seq_"+str(MAX)+".txt")
+    
+# Plot
+# plt.figure(figsize=(8, 6))
+# plt.scatter(lambda_2_values, convergence_times, color="blue", alpha=0.7)
+# plt.xlabel("Second Largest Eigenvalue (λ2)")
+# plt.ylabel("Convergence Time (Cycles)")
+# plt.title("Convergence Time vs. λ2")
+# plt.grid(True)
+# plt.show()

simu-opti/transition.py

@@ -0,0 +1,24 @@
+import numpy as np
+
+def converge(v,p):
+    i = 0
+    while (np.max(np.abs(v - np.mean(v))) > 1):
+        v = np.matmul(p, v)
+        i += 1
+    return i
+
+def create_vector(max):
+    np.random.seed(0)  # For reproducibility
+    vector = np.random.randint(0, 1000, size=(max, 1))  # Random integers
+    while np.any(np.abs(np.diff(vector.flatten())) < 100):  # Ensure the difference is at least 100
+        vector = np.random.randint(0, 1000, size=(max, 1))
+    return vector
+
+v = create_vector(3)
+p = np.array([
+    [0.500000, 0.250000, 0.250000],
+    [0.250000, 0.375000, 0.375000],
+    [0.000000, 0.250000, 0.750000]
+])
+
+print(converge(v,p))