1 points par GN⁺ 2024-07-05 | 1 commentaires | Partager sur WhatsApp
  • La multiplication de matrices de NumPy s’appuie sur des bibliothèques BLAS externes, mais cette implémentation vise à porter les performances mono-thread et multi-thread au niveau de BLAS avec seulement du C pur et FMA3·AVX2
  • Le cœur des performances consiste à découper $C$ en petits blocs et à utiliser un microkernel 16×6 qui répète des rank-1 updates dans les registres YMM afin de réduire les accès mémoire
  • Pour les matrices de taille arbitraire, la gestion des bords peut vite devenir un goulot d’étranglement ; l’implémentation combine donc des stores masqués et des tampons avec padding à 0 pour éviter la baisse de performance des loads masqués
  • La réutilisation du cache est assurée par un blocking k_c, m_c, n_c, et les performances maximales réelles dépendent fortement du tuning du nombre de threads, de la taille du kernel et de la taille des tuiles
  • AVX-512 a été exclu afin de prendre en charge un éventail plus large de CPU ; sur les CPU AVX-512, BLAS peut donc être plus rapide, et la comparaison avec OpenBLAS a elle aussi été réalisée avec AVX-512 désactivé

Objectif de l’implémentation et base de comparaison

  • Le code de l’implémentation est publié dans sgemm.c et optimise la multiplication de matrices FP32 multi-thread sur les processeurs récents
  • NumPy dépend de bibliothèques BLAS externes pour les opérations d’algèbre linéaire comme la multiplication de matrices
    • Exemples : Intel MKL, Accelerate, BLIS, GotoBLAS, OpenBLAS
    • OpenBLAS, GotoBLAS et BLIS sont écrits en C/FORTRAN/Assembly et incluent des implémentations de multiplication de matrices optimisées manuellement pour chaque microarchitecture CPU
  • L’objectif est une implémentation de multiplication de matrices écrite en C pur, sans assembleur bas niveau, tout en respectant les conditions suivantes
    • Fonctionner avec des matrices de tailles arbitraires
    • S’exécuter sur les processeurs x86-64 récents
    • Rivaliser avec les bibliothèques BLAS existantes
    • Rester simple et facile à étendre
  • Les références sont Fast Multidimensional Matrix Multiplication on CPU from Scratch de Simon Boehm, Matrix Multiplication de Sergey Slotin, Can you multiply a matrix? de Geohot, ainsi que des articles sur GotoBLAS et BLIS

Conditions de benchmark et calcul des FLOPS

  • L’environnement de test est AMD Ryzen 7 9700X, 32 Go DDR5 6000 MHz CL36, OpenBLAS 0.3.26, GCC 13.3, Ubuntu 24.04.1 LTS
  • Les flags de compilation utilisés sont -O3 -march=native -mno-avx512f -fopenmp
  • Pour une comparaison équitable, il faut définir le TARGET approprié lors de l’installation d’OpenBLAS et désactiver les instructions AVX-512
    • Les processeurs Zen4/5 se compilent avec make TARGET=ZEN
    • Sinon, OpenBLAS utilise les instructions AVX-512 par défaut
  • La multiplication de matrices FP32 d’OpenBLAS est exécutée via l’API cblas_sgemm
  • Le benchmark porte sur des matrices carrées
    • Évaluation de m=n=k=200 à m=n=k=10000, par pas de 200
    • La multiplication de matrices est répétée n_iter fois, et le temps d’exécution médian est utilisé pour mesurer les performances
  • Multiplier une matrice $A$ de taille $M \times K$ par une matrice $B$ de taille $K \times N$ représente un total de $2MNK$ FLOP
    • Les performances sont calculées avec FLOPS=(2*m*n*k)/exec_time

Limites théoriques et base SIMD

  • Les CPU x86-64 récents traitent plusieurs données en parallèle grâce aux extensions SIMD
  • Les principales instructions sont AVX2 et FMA
    • Toutes deux utilisent des registres YMM de 256 bits
    • Chaque registre YMM peut contenir 8 floats de 32 bits
  • L’instruction FMA VFMADD231PS effectue une opération packed single de la forme YMM1 = YMM2 * YMM3 + YMM1
  • Sur Ryzen 9700X, le débit du fused multiply-add est de 0,5 cycle/instruction, soit 2 instructions par cycle
  • En théorie, le Ryzen 9700X peut effectuer 32 FLOP par cycle sur un seul cœur
    • Le calcul est 8 floats × 2(add+mul) × 2(1/TP)
    • En supposant une fréquence soutenue de 4,7 GHz sur 8 cœurs, le pic théorique multi-thread est estimé à 1203 FLOPS

Implémentation de base et microkernel

  • Les matrices sont stockées en ordre column-major
    • A[row][col] est accessible, depuis un pointeur C, via ptr[col*M + row]
  • L’implémentation la plus simple parcourt toutes les lignes et colonnes de $C$ et calcule, pour chaque élément, le produit scalaire d’une ligne de $A$ et d’une colonne de $B$
  • Le cœur d’une implémentation haute performance est le microkernel, qui découpe $C$ en sous-matrices $m_R \times n_R$ et calcule efficacement chacune d’elles
  • Le kernel initialise $\bar{C}$ à zéro dans les registres, puis itère le long de la dimension $K$
    • Il charge dans les registres le vecteur colonne de $\bar{A}$ et le vecteur ligne de $\bar{B}$
    • Il calcule le produit extérieur de ces deux vecteurs et l’ajoute aux accumulateurs $\bar{C}$
    • Chaque étape est une rank-1 update
  • Par rapport aux $2K m_R n_R$ accès mémoire de l’approche naïve, cette méthode réduit à $(m_R+n_R)K$ le nombre d’éléments chargés dans les registres
  • Les CPU AVX disposent de 16 registres YMM, la taille du kernel doit donc respecter la contrainte suivante
    • $(m_R/8) \cdot n_R + m_R/8 + 1 \le 16$
    • $m_R$ doit être un multiple de 8
  • En théorie, plus $m_R$ et $n_R$ sont grands et proches l’un de l’autre, plus la réduction des accès mémoire est importante, mais en pratique, sur Ryzen 9700X, le kernel 16×6 donne les meilleures performances
  • L’implémentation utilise les intrinsics de immintrin.h
    • __m256 est un type vectoriel 256 bits représentant le contenu d’un registre YMM
    • _mm256_loadu_ps charge un vecteur colonne de A
    • _mm256_broadcast_ss diffuse une valeur scalaire de B dans un vecteur de 8 floats
    • _mm256_fmadd_ps met à jour les accumulateurs
    • _mm256_storeu_ps stocke le résultat en mémoire
  • L’assembleur généré contient des instructions SIMD FMA comme vfmadd231ps et vbroadcastss

Padding pour matrices de taille arbitraire

  • Le kernel 16×6 de base fonctionne directement lorsque $M$ et $N$ sont respectivement des multiples de 16 et de 6
  • Sur les zones de bord, si le nombre de colonnes $n$ est inférieur à 6, la boucle de stockage ne s’exécute que jusqu’à j < n
  • Lorsque le nombre de lignes $m$ est inférieur à 16, _mm256_storeu_ps stocke 8 éléments à la fois ; un store masqué est donc nécessaire
    • _mm256_maskstore_ps n’écrit en mémoire que les éléments dont le bit de masque est activé
    • Le masque est généré selon le nombre de lignes recouvrées $m$
  • Traiter aussi les loads de bord avec _mm256_maskload_ps peut fortement dégrader les performances du kernel
    • Les instructions supplémentaires de calcul du masque ajoutent un surcoût
    • Comme $n$ n’est pas une constante connue à la compilation, le compilateur a du mal à dérouler efficacement la boucle
  • À la place, si $m \neq m_R$, $\bar{A}$ est copié dans un tampon et complété par des zéros ; si $n \neq n_R$, $\bar{B}$ est également copié dans un tampon et rempli de zéros
  • L’implémentation correspondante se trouve dans matmul_pad.h

Blocking de cache et réutilisation des données

  • Entre les registres et la DRAM se trouve la hiérarchie de cache CPU ; les CPU de bureau récents utilisent généralement des caches L1, L2 et L3
  • Le cache est plus rapide que la DRAM mais sa capacité est limitée, il est donc impossible de placer l’intégralité de $A$, $B$ et $C$ en cache
  • La méthode consistant à découper les matrices en petits blocs, à les charger dans le cache et à réutiliser les mêmes données pour plusieurs rank-1 updates s’appelle blocking de cache ou tiling
  • Le blocking de cache mono-thread prend la forme d’une boucle quintuple similaire à la structure de BLIS
    • La boucle la plus externe parcourt la dimension $N$ pour créer les blocs $C_j$ et $B_j$
    • La boucle suivante parcourt la dimension $K$ pour créer les blocs $A_j$ et $B_p$
    • $B_p$ est packé pour devenir $\tilde{B}_p$ et, si nécessaire, complété par des zéros afin de favoriser sa réutilisation dans le cache L3
    • La boucle suivante parcourt la dimension $M$ pour créer les blocs $C_i$ et $A_j$, et $A_j$ est packé pour devenir $\tilde{A}_j$
    • Les deux dernières boucles découpent les blocs de cache en panneaux $m_R \times k_c$ et $k_c \times n_R$, puis les transmettent au kernel
  • Les matrices packées $\tilde{A}_j$ et $\tilde{B}_p$ sont stockées différemment
    • Les panneaux internes de $\tilde{A}_j$ sont stockés en column-major
    • Les panneaux internes de $\tilde{B}_p$ sont stockés en row-major
  • Les paramètres de blocking de cache doivent être ajustés aux tailles de cache de chaque modèle de CPU
    • $k_c \times n_c$ est un point de départ pour remplir le cache L3
    • $m_c \times k_c$ est un point de départ pour remplir le cache L2
    • $k_c \times n_R$ est un point de départ pour remplir le cache L1
  • En pratique, des valeurs supérieures aux valeurs théoriques donnent souvent de meilleures performances, et comme le CPU gère automatiquement le placement en cache, il faut concevoir les boucles et les schémas d’accès au niveau algorithmique
  • L’implémentation se trouve dans matmul_cache.h

Micro-optimisations du kernel

  • Au lieu de définir les accumulateurs sous forme de tableau, comme __m256 C_buffer[6][2], les variables d’accumulateurs sont explicitement déroulées et déclarées séparément
  • Cette approche aide GCC à mieux optimiser le code et à éviter le register spilling
  • Le calcul du masque est également modifié pour utiliser des instructions vectorielles
    • Un tableau statique mask[32] est utilisé avec _mm256_cvtepi8_epi32 et _mm_loadu_si64
  • Cette implémentation se trouve dans matmul_micro.h

Stratégie multi-threading

  • La parallélisation vise à la fois les opérations arithmétiques et le packing
  • Les 5e, 4e et 3e boucles en partant du microkernel répètent des blocs de la taille des blocs de cache
    • Pour occuper tous les threads, le nombre d’itérations doit être au moins égal au nombre de threads
    • Les dimensions des matrices d’entrée doivent être approximativement au moins égales à nombre de threads × taille de bloc de cache
  • Sur Ryzen 9700X, les tailles de blocs de cache qui donnent de bonnes performances en mono-thread sont $n_c=1535$, $m_c=1024$
    • Pour exploiter les 8 cœurs, il faut au minimum une dimension de taille $\max(m_c,n_c) \times 8 = 1535 \times 8 = 12280$
  • À l’inverse, les deux dernières boucles répètent de petits blocs $m_R$, $n_R$ et conviennent donc à la parallélisation
    • En général, $m_R$ et $n_R$ sont inférieurs à 20
    • Choisir $m_c$ et $n_c$ comme multiples du nombre de cœurs permet de répartir équitablement le travail
  • Sur Ryzen 9700X, la meilleure performance est obtenue en parallélisant ensemble les deux boucles internes avec #pragma omp parallel for collapse(2) num_threads(NTHREADS)
  • Sur les processeurs avec beaucoup de cœurs, en particulier au-delà de 16 cœurs, on peut envisager le parallélisme imbriqué et la parallélisation de 2 à 3 boucles
  • Le packing de $\tilde{A}$ et $\tilde{B}$ est également parallélisé avec OpenMP
    • pack_blockA est parallélisé en parcourant mc par pas de MR
    • pack_blockB est parallélisé en parcourant nc par pas de NR
  • Les paramètres qui donnent de bonnes performances sur Ryzen 9700X dans l’implémentation multi-thread sont les suivants
    • $m_c = m_R \times \text{nombre de threads} \times 5$
    • $n_c = n_R \times \text{nombre de threads} \times 50$
  • L’implémentation multi-thread finale se trouve dans matmul_parallel.h

1 commentaires

 
GN⁺ 2024-07-05
Commentaires sur Hacker News
  • Si l’idée de cet article est simplement qu’il reste généralement de la marge de performance, alors il sous-estime plutôt l’ampleur possible des améliorations. Et ce, même si l’effort investi dans les bibliothèques de multiplication de matrices est bien plus important que dans la plupart des logiciels.
    À moins d’avoir déjà un code fortement optimisé, il est courant d’obtenir, sans énorme effort, des gains de 10 à 1000 fois, voire plus par rapport au code existant. En gros, par ordre d’importance, le plus crucial est de savoir si le choix de l’algorithme est approprié et si l’on peut supprimer purement et simplement le travail à effectuer ; réduire les opérations coûteuses comme les allers-retours vers le noyau ou les malloc compte aussi beaucoup.
    Pour la vectorisation, les intrinsics vectoriels explicites sont utiles, mais le simple fait de réorganiser les données, par exemple en passant d’un tableau de structures à une structure de tableaux, produit souvent le même code machine. L’efficacité du cache est également importante, et dans le code parallèle cela devient plus complexe lorsque les données ne sont pas bien isolées par thread, comme avec le faux partage. Enfin, on peut aussi faire des optimisations spécifiques au matériel, via des intrinsics ou de l’assembleur écrit à la main.

    • Il ne faut pas non plus négliger l’impact du réseau. Une fois, j’ai découvert qu’une requête distribuée rapatriait environ un million de lignes par le réseau, puis faisait une jointure pour n’en garder que 5 à 10 ; corriger cela a permis un gain de performance de plusieurs centaines de fois.
      En modifiant la requête pour que la jointure se fasse sur le serveur distant et que seules 5 à 10 lignes soient envoyées sur le réseau, c’est devenu immédiatement plus rapide. Il y a toujours un surcoût fixe et de la latence, mais si l’on envoie beaucoup plus de données que nécessaire sur une connexion réseau, les performances finissent par s’effondrer. L’article « It's the latency, stupid », sur l’impact de la latence, vaut aussi la lecture : http://www.stuartcheshire.org/rants/latency.html
      Dans l’ensemble, je suis d’accord avec les considérations ci-dessus et leur ordre approximatif.
    • « Le choix de l’algorithme est-il approprié ? » est en pratique devenu une forme de cargo cult. Un algorithme « plus rapide » peut avoir des facteurs constants horribles, et il arrive souvent que l’approche qui fait davantage de travail soit en réalité plus performante.
      Beaucoup d’entretiens sont devenus des quiz de mémorisation d’algorithmes obscurs façon « Google fait comme ça », au lieu d’évaluer la capacité à raisonner sur les causes de lenteur d’une implémentation, à la benchmarker et à la corriger.
  • Les schémas de code courants ne sont pas assez spécialisés pour le matériel et laissent beaucoup de performance sur la table. Cet article en est un exemple intéressant ; une autre démonstration classique est « There's plenty of room at the top ».
    https://www.science.org/doi/10.1126/science.aam9744

  • Pour comprendre ce sujet, les articles du dépôt BLIS sont ce qui se rapproche le plus d’une référence. Je ne vois pas pourquoi on penserait qu’un BLAS optimisé ne donne pas de bonnes performances ; pour des matrices suffisamment grandes, il faut s’attendre à plus de 90 % du pic du CPU.
    La dernière fois que j’ai regardé, OpenBLAS en série était globalement comparable à MKL, et BLAS implémente GEMM comme brique de base d’algèbre linéaire, pas matmul. Je ne comprends pas non plus pourquoi on utiliserait numpy au lieu d’un framework de benchmark, et sur Zen il faudrait à mon avis comparer avec le BLAS d’AMD, c’est-à-dire l’implémentation basée sur BLIS. BLIS avait autrefois une meilleure histoire côté parallélisation qu’OpenBLAS, et AMD BLIS dispose aussi d’un basculement d’implémentation pour les « petites » dimensions ; je ne sais pas si OpenBLAS l’a aujourd’hui.
    Les intrinsics SIMD ne sont pas indispensables pour vectoriser un micro-kernel ; un bon compilateur C peut vectoriser complètement et dérouler les boucles. Le micro-kernel en C pur de BLIS atteint plus de 80 % des performances de l’implémentation optimisée à la main pour Haswell, avec des tailles de blocs appropriées. La différence vient probablement du préchargement, mais je ne la comprends pas précisément.

    • Les intrinsics SIMD et le déroulage manuel des boucles sont clairement nécessaires. C’est pour cela que toutes les bibliothèques BLAS vectorisent et déroulent les boucles à la main.
      Même les compilateurs modernes ne réussissent pas l’auto-vectorisation et le déroulage de boucles correctement dans 100 % des cas.
  • L’article et l’implémentation ont l’air bons, mais je me demande quel est le « secret ». OpenBLAS est optimisé en assembleur+C depuis des décennies précisément pour ce type de problème ; comment peut-il être battu ?
    L’article détaille le cache, etc., mais je me demande si BLAS n’exploite pas ce genre de choses, ou si c’est simplement mieux adapté à un processeur particulier.

    • OpenBLAS n’est pas optimisé à ce point pour certaines architectures récentes. Les matrices n’étaient pas non plus si grandes, et numpy a le surcoût de cffi.
      L’écart de performance était bien plus marqué sur le débit de pointe que sur le débit moyen, or il y a très peu d’applications où le pic compte vraiment. Le code de benchmark affiché semble faire passer numpy par l’allocateur Python, alors que l’implémentation C ne passe pas par l’allocateur ; c’est donc le premier endroit où chercher une erreur ou une incohérence de microbenchmark. Beaucoup de routines numpy prennent en charge les opérations en place, donc il faudrait explicitement voir un benchmark des versions en place des deux côtés.
      numpy effectue aussi des vérifications de bornes et de la gestion d’erreurs indépendantes de l’implémentation sous-jacente, ce qui explique pourquoi il peut être très lent sur de petites matrices, parfois même plus lent que des listes Python brutes. Ajouter quelques milliers de cycles de pur overhead rend vite les choses difficiles à accélérer.
      Cette implémentation suit une approche assez méthodique pour saturer les caches concernés et, en un sens, c’est évident, mais des améliorations d’ingénierie claires comme celle-ci méritent d’être mises en avant dans ce type de discussion. OpenBLAS a aussi mobilisé beaucoup de monde, mais il est peu probable qu’ils aient pensé à tout. Pour l’expliquer correctement, il faudrait une analyse approfondie du code des deux côtés.
    • Battre OpenBLAS n’a rien de surprenant ni d’inédit. Par exemple, la bibliothèque d’algèbre linéaire Mir du langage D l’a déjà fait il y a quelques années [1].
      Pour des implémentations C++ et C, voir les approches par métaprogrammation [2], [3]. Ce qui est vraiment surprenant, c’est que beaucoup de langages modernes comme Matlab, Julia ou Mojo dépendent encore d’OpenBLAS, même s’ils ont bien sûr chacun leurs raisons.
      [1] Numeric age for D: Mir GLAS is faster than OpenBLAS and Eigen (2016):
      http://blog.mir.dlang.io/glas/benchmark/openblas/2016/09/23/...
      [2] Vastly outperforming LAPACK with C++ metaprogramming (2018):
      https://wordsandbuttons.online/vastly_outperforming_lapack_w...
      [3] Outperforming LAPACK with C metaprogramming (2018):
      https://wordsandbuttons.online/outperforming_lapack_with_c_m...
    • -march=native peut donner un avantage, car il compile pour le modèle exact de CPU. numpy a de fortes chances d’avoir été compilé pour une cible x86-64 plus générique et plus ancienne.
      Sur les CPU Ryzen, -march=native utilise probablement v4, tandis que numpy cible sans doute v1 ou v2.
      https://en.wikipedia.org/wiki/X86-64#Microarchitecture_level...
    • numpy 2.0 intègre Google Highway afin de mieux utiliser SIMD sur plusieurs microarchitectures, donc la comparaison devrait devenir plus favorable à numpy.
  • L’article est bon, et c’est aussi très bien d’avoir rendu les benchmarks faciles à reproduire. Sur mon Xeon W-2245 16 cœurs à 3,90 GHz, matmul.c effectue une multiplication de matrices 8192x8192 en 1,41 s avec gcc -O3 et en 1,47 s avec clang -O2, tandis que NumPy prend 1,07 s.
    Avec un kernel AVX-512, je pense que ce serait beaucoup plus rapide. Une autre raison possible des performances décevantes est OpenMP ; d’après mon expérience, gérer explicitement un pool de threads avec pthreads permet de réduire l’overhead. Il vaudrait aussi mieux utiliser sysconf(_SC_NPROCESSORS_ONLN) plutôt que de coder en dur le nombre de CPU.

  • Il n’y a pas de raison d’imposer des charges différentes, avec Python d’un côté et C de l’autre. On aurait pu écrire les deux en C, l’un appelant une bibliothèque BLAS et l’autre cette implémentation, pour faire une comparaison à armes égales.

    • Ici, la comparaison avec Python est pertinente. Aujourd’hui, la manière la plus populaire d’effectuer ce type de calcul est d’utiliser Python avec numpy.
      L’overhead n’est pas énorme, mais comme dit ailleurs dans ce fil, il est important de l’appeler correctement. Opposer du code numpy naïf à du code C optimisé n’est clairement pas une comparaison équitable.
  • Ce n’est pas un chemin critique, mais l’inefficacité de la génération du masque, c’est-à-dire l’utilisation de bit_mask, me gêne. Des méthodes plus efficaces consisteraient à créer un tableau constant global de la forme {-1,-1,...,0,0,...} et à charger à partir des offsets d’éléments 16-m, 8-m, ou à comparer un vecteur constant {0,1,2,3,4,...} avec m et m-8 diffusés par broadcast
    Cela dit, cela ne concerne qu’une seule colonne de la matrice, et la boucle maskload/maskstore qui suit prend beaucoup plus de temps, donc c’est vraiment un détail. En particulier, les stores restent lents même sur Zen 4[1], et les instructions AVX-512 sont 6 fois plus rapides alors que la seule différence est qu’elles prennent le masque depuis un registre de masque. Clang vectorise de toute façon automatiquement les décalages, donc ce sera probablement seulement 2 à 3 fois plus lent que ma suggestion
    [1]: https://uops.info/table.html?search=vmaskmovps&cb_lat=on&cb_...

    • Je suis l’auteur. C’est vraiment la première fois que j’optimise du code C et que j’utilise des intrinsics, donc je ne suis pas expert du domaine, mais j’aimerais en apprendre davantage
      J’apprécie vraiment les retours qui apportent un nouveau point de vue. Si je me souviens bien, « créer un tableau global de constantes et le charger » était un peu plus lent que le décalage du masque de bits quand je l’avais testé, mais je vais retester pour en être sûr. La méthode consistant à « comparer un vecteur constant {0, 1, 2, 3, 4, ...} avec m et m-8 diffusés par broadcast » est une bonne idée, je vais l’essayer
    • Pour créer le tableau global de constantes, on peut mettre les éléments en int8_t et, lors du chargement, faire une extension de signe des octets en int32_t. La combinaison _mm_loadu_si64 / _mm256_cvtepi8_epi32 devrait être compilée en une seule instruction vpmovsxbd avec opérande mémoire
      Ainsi, avec un alignement correct via alignas(32), tout le tableau de constantes tient dans une seule ligne de cache. Dans le cas d’usage de l’article, il faut deux masques, donc la deuxième instruction vpmovsxbd sera clairement un hit dans le cache L1D, ce qui convient bien
  • Et tinyBLAS de jart ?
    https://hacks.mozilla.org/2024/04/llamafiles-progress-four-m...
    Et https://justine.lol/matmul/

    • J’ai beaucoup discuté avec Justine hier, et sur cette station de travail, cette implémentation semble être au moins 2 fois plus rapide que tinyBLAS. Toute la discussion est sur le Discord Mozilla AI : https://discord.com/invite/NSnjHmT5xY
  • À part les benchmarks, quelle est la raison de multithreader la multiplication matricielle elle-même ? En pratique, ne vaudrait-il pas mieux utiliser plusieurs threads au niveau de l’algorithme qui utilise la multiplication ?

    • En HPC, c’est effectivement généralement ce qu’on fait. Cela dit, remplacer simplement par un BLAS parallèle peut facilement aider certains types de code R
      Mais dans le code HPC, GEMM n’est généralement pas le goulot d’étranglement
  • Je n’ai encore fait que survoler, mais cet article contient beaucoup de détails et d’explications. Il semble assez bien expliquer comment une multiplication matricielle rapide est implémentée en tenant compte des considérations d’architecture, donc je l’ai ajouté à ma liste de lecture