Élever une matrice à une puissance, en Python

Exercice d'entraînement · niveau 3 (difficile) · mathématiques appliquées (ECG 1re année), chapitre 2 — Calcul matriciel et systèmes linéaires · Puissances et polynômes annulateurs

Énoncé

Écrire une fonction Python puissance(A, n) qui calcule $A^{n}$ pour une matrice carrée A (tableau numpy) et un entier $n \geq 0$ , en n'utilisant que np.dot, np.eye et une boucle. Combien de produits de matrices effectue-t-elle ? Proposer ensuite une version qui en fait beaucoup moins en utilisant $A^{2 k} = (A^{k})^{2}$ .

Corrigé

Version directe.

import numpy as np

def puissance(A, n):
    P = np.eye(len(A))          # A^0 = I
    for _ in range(n):
        P = np.dot(P, A)        # un produit par tour
    return P

La boucle tourne $n$ fois et fait un produit par tour : $n$ produits de matrices. Pour des matrices $p \times p$ , chaque produit coûte de l'ordre de $p^{3}$ multiplications de nombres, soit environ $n p^{3}$ en tout.

Version rapide. On exploite $A^{2 k} = (A^{k})^{2}$ et $A^{2 k + 1} = A \times (A^{k})^{2}$ :

def puissance_rapide(A, n):
    if n == 0:
        return np.eye(len(A))
    B = puissance_rapide(A, n // 2)
    C = np.dot(B, B)
    if n % 2 == 1:
        C = np.dot(A, C)
    return C

Le compte. À chaque appel, $n$ est divisé par $2$ : il y a donc de l'ordre de $lo g_{2} n$ appels, et chacun fait un ou deux produits. Le nombre de produits est au plus $2 lo g_{2} n$ , contre $n$ pour la première version. Pour $n = 1000$ : environ $20$ produits au lieu de $1000$ .

Les autres exercices de ce chapitre Le cours du chapitre

Un blocage sur cet exercice ? Le tuteur d'Adloun guide par questions, sans donner la réponse.