#------------------------TD------------------------#


#####################[Code 1]##########################

import numpy as np

def f(M, s, N):
    return (M + s, M * s, M ** 2, M + N, M - N, M * N)

M = np.array([[1, 2], [3, 4]])
N = np.array([[5, 6], [7, 8]])
s = 10

plus_scal, fois_scal, carre, plus_mat, moins_mat, produit = f(M, s, N)

print("M + s:\n", plus_scal)
print("M * s:\n", fois_scal)
print("M ** 2:\n", carre)
print("M + N:\n", plus_mat)
print("M - N:\n", moins_mat)
print("M * N:\n", produit)



#####################[Code 2]##########################

import numpy as np

A = np.array([[1, 2], [3, 4]])
B = np.array([[5, 6], [7, 8]])

print("A^T:\n", A.T, "\n")
print("A+A:\n", A + A, "\n")
print("AA^T:\n", A @ A.T, "\n")
print("A^TA:\n", A.T @ A, "\n")
print("AB:\n", A @ B, "\n")

def f(A, B, mylambda):
    I = np.identity(B.shape[0])
    return A @ (B - mylambda * I)

print("A(B - lambda I) avec lambda =2:\n", f(A, B, 2))



#####################[Code 3]##########################

import numpy as np

def cosinus_sim(u, v):
    produit = np.dot(u, v)
    norm_u = np.linalg.norm(u)
    norm_v = np.linalg.norm(v)

    if norm_u == 0 or norm_v == 0:
        return 0.0
    return produit / (norm_u * norm_v)

u1 = np.array([1, 2, 3])
v1 = np.array([4, 5, 6])
print(f"cos({u1}, {v1}) = {cosinus_sim(u1, v1):.4f}")


u4 = np.array([1, 1])
v4 = np.array([-1, -1])
print(f"cos({u4}, {v4}) = {cosinus_sim(u4, v4):.4f}")

u5 = np.array([0, 0, 0])
v5 = np.array([1, 2, 3])
print(f"cos(zero, v) = {cosinus_sim(u5, v5)}")



#####################[Code 4]##########################

import numpy as np

def normale_triangle(p1, p2, p3):
    u = p2 - p1
    v = p3 - p1
    n = np.cross(u, v)
    return n / np.linalg.norm(n)

p1 = np.array([0, 0, 0])
p2 = np.array([1, 0, 0])
p3 = np.array([0, 1, 0])
print("Normale :", normale_triangle(p1, p2, p3))

p1 = np.array([1, 0, 0])
p2 = np.array([0, 1, 0])
p3 = np.array([0, 0, 1])
print("Normale :", normale_triangle(p1, p2, p3))



#####################[Code 5]##########################

import numpy as np

def normes_carrees(A):
    return np.sum(A**2, axis=1)

A = np.array([[1, 2], [3, 4], [5, 6]])
print(normes_carrees(A))



#####################[Code 6]##########################

import numpy as np

def distances_point_centres(x, A):
    diff = x - A   #<-- broadcasting ici!
    return np.sum(diff**2, axis=1)

A = np.array([[0, 0], [1, 1], [2, 2]])
x = np.array([1, 1])
print( distances_point_centres(x, A) )



#####################[Code 7]##########################

import numpy as np

def g(A, L):
    indices = np.where(np.abs(A) >= L)[0]

    return indices[0] if indices.size > 0 else None

print("Exemple 1:", g(np.array([0, 1, 2, 3, 4, 5]), 3))
print("Exemple 2:", g(np.array([0, 1, 2, 1, 0, -1]), 3))
print("Exemple 3:", g(np.array([0, -1, -2, -3, -4]), 3))
print("Exemple 4:", g(np.array([5, 4, 3, 2, 1]), 3))
print("Exemple 5:", g(np.array([1, -3, 2, 4]), 3))



#####################[Code 8]##########################

import numpy as np

def findClosest(A, z):
    return A[np.argmin(np.abs(z - A))]

print( findClosest( np.array([2,1,5,1,22]), 1) )



#####################[Code 9]##########################

import numpy as np

def analyse(A):
    return np.max(A, axis=1), np.any(A, axis=1)

A = np.array([[False, True, False],
              [False, False, False],
              [True, False, True]])
print(analyse(A))



#####################[Code 10]##########################

import numpy as np

def auto_corr_circ(a):
    n = len(a)
    c = np.zeros(n)
    for j in range(n):
        c[j] = np.dot(a, np.roll(a, j))
    return c

print( auto_corr_circ(np.array([1,2,3])) )



#####################[Code 11]##########################

import numpy as np

def centre_points(points):
    return np.mean(points, axis=0)

points = np.array([[1, 2], [3, 4], [5, 6]])
print(centre_points(points))



#####################[Code 12]##########################

import numpy as np

def proche(a, b, eps=1e-4):
    return np.allclose(a, b, atol=eps, rtol=0)

a = np.array([1.00001, 2.00001, 3.00001])
b = np.array([1.0, 2.0, 3.0])
print( proche(a, b, eps=1e-4) )
print( proche(a, b, eps=1e-6) )



#####################[Code 13]##########################

import numpy as np

def evalue_polynome(x):
    coeffs = [3, -2, 0, 5]
    return np.polyval(coeffs, x)

print(evalue_polynome(np.array( [1,2,3] )))



#####################[Code 14]##########################

import numpy as np

def racines():
    return np.roots([1, 0, -5, 0, 4])

print(racines())



#####################[Code 15]##########################

import numpy as np

def insertion_triee(arr, x):
    idx = np.searchsorted(arr, x, side='right')
    return np.insert(arr, idx, x)

# Exemples
print(insertion_triee(np.array([1, 3, 5]), 4))    # [1, 3, 4, 5]
print(insertion_triee(np.array([1, 3, 3, 5]), 3)) # [1, 3, 3, 3, 5]
print(insertion_triee(np.array([1, 2, 3]), 0))    # [0, 1, 2, 3]
print(insertion_triee(np.array([1, 2, 3]), 4))    # [1, 2, 3, 4]



#####################[Code 16]##########################

import numpy as np

def mode(arr):
    valeurs, counts = np.unique(arr, return_counts=True)
    return valeurs[np.argmax(counts)]

print(mode(np.array([1, 2, 1, 3, 1, 2])))
print(mode(np.array([1, 1, 2, 2, 3])))
print(mode(np.array([5, 5, 5, 2, 2, 2, 3])))
print(mode(np.array(['a', 'b', 'a', 'c', 'a', 'b'])))



#####################[Code 17]##########################

import numpy as np

def random_operations(N=100, M=50, k=5):
    A = np.random.choice([-1, 1], size=(N, M))

    sommes = np.sum(A, axis=1)

    indices = np.random.choice(N, k, replace=False)

    return A, sommes, indices

A, sommes, indices = random_operations(10, 5, 3)
print("Matrice A:\n", A)
print("Sommes par ligne:", sommes)
print("Indices selectionnes:", indices)



#####################[Code 18]##########################

import numpy as np

def f(A):
    positions = np.cumsum(A, axis=1)

    positions_finales = positions[:, -1]

    max_abs = np.max(np.abs(positions), axis=1)

    return positions_finales, max_abs

A1 = np.array([ [1, 1, 1, 1], [-1, 1, -1, 1]])
finales, max_abs = f(A1)
print("A1:\n", A1)
print("Positions finales:", finales)
print("Max valeur absolue:", max_abs)

np.random.seed(42)
A2 = np.random.choice([-1, 1], size=(5, 10))
finales, max_abs = f(A2)
print("A2:\n", A2)
print("Positions finales:", finales)
print("Max valeur absolue:", max_abs)



#####################[Code 19]##########################

import numpy as np

def mult_blocs(A, B, bs=64):
    n = A.shape[0]
    C = np.zeros((n, n))
    for i in range(0, n, bs):
        for j in range(0, n, bs):
            for k in range(0, n, bs):
                i_end = min(i+bs, n)
                j_end = min(j+bs, n)
                k_end = min(k+bs, n)
                C[i:i_end, j:j_end] += A[i:i_end, k:k_end] @ B[k:k_end, j:j_end]
    return C


A = np.random.random( (3,3) )
B = np.random.random( (3,3) )
print(A@B)
print(mult_blocs(A,B))



#####################[Code 20]##########################

import numpy as np

def randomwalk(L=10, n_sim=10000, max_steps=100000):
    A = np.random.choice([-1, 1], size=(n_sim, max_steps))
    positions = np.cumsum(A, axis=1)

    temps = []
    for i in range(n_sim):
        indices = np.where( np.abs(positions[i]) == L )[0]

        if len(indices) > 0:
            temps.append(indices[0] + 1)
        else:
            temps.append(max_steps)

    return sum(temps) / len(temps)

print( randomwalk (3, 10, 12) )



#####################[Code 21]##########################

import numpy as np

def kmeans(X, k, max_iter=100):
    n, d = X.shape

    idx = np.random.choice(n, k, replace=False)
    centres = X[idx].copy()

    for iteration in range(max_iter):
        labels = np.zeros(n)
        for i in range(n):
            distances = np.zeros(k)
            for j in range(k):
                distances[j] = np.sum( (X[i] - centres[j]) ** 2)
            labels[i] = np.argmin(distances)

        nouveaux_centres = np.zeros((k, d))
        for j in range(k):
            indices_j = np.where(labels == j)[0]
            if len(indices_j) > 0:
                nouveaux_centres[j] = np.mean(X[indices_j], axis=0)
            else:
                nouveaux_centres[j] = centres[j]

        if np.allclose(centres, nouveaux_centres, atol=1e-4):
            print(f"Convergence apres {iteration + 1} iterations")
            break
        centres = nouveaux_centres

    return labels, centres


# ============================================================
# EXEMPLE 1 :     (4 points, 2 clusters)
# ============================================================
print("=" * 50)
print("EXEMPLE 1 : 4 points, 2 clusters")
print("=" * 50)

X1 = np.array([
    [1, 1],   # point 0
    [1, 2],   # point 1
    [8, 8],   # point 2
    [9, 8]    # point 3
])

print("Points:")
print(X1)
print()

labels1, centres1 = kmeans(X1, k=2, max_iter=10)

print(f"\nLabels: {labels1}")
print(f"Centre cluster 0: {centres1[0]}")
print(f"Centre cluster 1: {centres1[1]}")


# ============================================================
# EXEMPLE 2 :  (2 clusters naturels)
# ============================================================
print("\n" + "=" * 50)
print("EXEMPLE 3 : Points sur un cercle (2 cotes)")
print("=" * 50)

np.random.seed(42)
theta1 = np.random.uniform(0, np.pi, 30)
theta2 = np.random.uniform(np.pi, 2*np.pi, 30)
r = 5

X3_haut = np.array([r * np.cos(theta1), r * np.sin(theta1)]).T
X3_bas = np.array([r * np.cos(theta2), r * np.sin(theta2)]).T
X3 = np.vstack([X3_haut, X3_bas])

labels3, centres3 = kmeans(X3, k=2, max_iter=20)

print(f"\nRepartition:")
print(f"  Cluster 0: {np.sum(labels3 == 0)} points")
print(f"  Cluster 1: {np.sum(labels3 == 1)} points")
print(f"\nCentres:")
print(f"  Centre 0: ({centres3[0][0]:.2f}, {centres3[0][1]:.2f})")
print(f"  Centre 1: ({centres3[1][0]:.2f}, {centres3[1][1]:.2f})")



# ============================================================
# Fonction   pour visualiser
# ============================================================
def plot_clusters(X, labels, centres, title="k-means Clustering"):
    import matplotlib.pyplot as plt

    colors = ['red', 'blue', 'green', 'orange', 'purple', 'brown', 'pink', 'gray']

    plt.figure(figsize=(8, 6))

    for j in range(centres.shape[0]):
        mask = labels == j
        plt.scatter(X[mask, 0], X[mask, 1],
                    c=colors[j % len(colors)],
                    label=f'Cluster {j}',
                    alpha=0.6, s=50)

    plt.scatter(centres[:, 0], centres[:, 1],
                c='black', marker='X', s=200,
                edgecolors='white', linewidth=2,
                label='Centres')

    plt.title(title)
    plt.legend()
    plt.axis('equal')
    plt.show()

plot_clusters(X3, labels3, centres3, "k-means sur deux demi-cercles")



#####################[Code 22]##########################

import numpy as np

def moyenne_cluster(X, labels, j):
    points = X[labels == j]
    if len(points) > 0:
        return np.mean(points, axis=0)
    else:
        return np.zeros(X.shape[1])

X = np.array([[1, 2], [3, 4], [5, 6], [7, 8]])
labels = np.array([0, 1, 0, 1])
print(moyenne_cluster(X, labels, 0))
print(moyenne_cluster(X, labels, 2))



#####################[Code 23]##########################

import numpy as np

def toeplitz(c, r):
    n = len(c)
    T = np.zeros((n, n))
    np.fill_diagonal(T, c[0])

    for k in range(1, n):
        T[0, k:] = r[k:]
        T[k, k:] = r[:n-k]

    for k in range(1, n):
        T[k:, 0] = c[k:]
        T[k:, k] = c[:n-k]

    return T

print(toeplitz(np.array([1,2,3]), np.array([1,4,5])))



#####################[Code 24]##########################

import numpy as np

def somme_sur_diagonale(A):
    return np.trace(A, offset=1)

print(somme_sur_diagonale(np.array([ [1,2,3], [4,5,6], [7,8,9] ])))



#####################[Code 25]##########################

def rang_colonne_plein(A):
    return np.linalg.matrix_rank(A) == A.shape[1]



#####################[Code 26]##########################

def somme_carres(m, n):
    I, J = np.indices((m, n))
    return I**2 + J**2



#####################[Code 27]##########################

import numpy as np

def randomwalk(L=10, n_sim=10000, max_steps=100000):
    A = np.random.choice([-1, 1], size=(n_sim, max_steps))
    positions = np.cumsum(A, axis=1)

    temps = []
    for i in range(n_sim):
        indices = np.where( np.abs(positions[i]) >= L )[0]

        if len(indices) > 0:
            temps.append(indices[0] + 1)
        else:
            temps.append(max_steps)

    return sum(temps) / len(temps)

print( randomwalk (3, 10, 12) )



#####################[Code 28]##########################

def kmeans(X, k, max_iter=100):
    n = len(X)
    idx = np.random.choice(n, k, replace=False)
    centres = X[idx]
    for _ in range(max_iter):
        diff = X[:, None, :] - centres[None, :, :]
        distances = np.sum(diff**2, axis=2)
        labels = np.argmin(distances, axis=1)
        nouveaux = np.array([
            X[labels == j].mean(axis=0) if (labels == j).any() else centres[j]
            for j in range(k)
        ])
        if np.allclose(centres, nouveaux):
            break
        centres = nouveaux
    return labels, centres



