Fondamentaux du Machine Learning 2.5 : Comportement générique en C avec les macrosDraft

2 octobre 2023

Le compromis avec la compréhension acquise en utilisant un langage de bas niveau comme C est le manque de fonctionnalités. Lorsqu’on écrit du code à des fins de démonstration, ce n’est pas vraiment un problème. Mais quand on doit écrire une application réelle, ça l’est.

L’une des fonctionnalités notables absentes du C est la généricité à la compilation.

La généricité à la compilation est un type de duplication de code qui dépend des contraintes suivantes :

Il existe plusieurs classes qui implémentent le même trait
Il existe du code qui ne dépend de ces classes qu’à travers le trait implémenté

Dans cet article, je veux vous montrer une façon de reproduire la généricité en C, en particulier par rapport à un trait d’opération dont dépendront toutes les opérations matricielles élément par élément. Nous allons faire cela avec la seule méthode intégrée de duplication de code en C : les macros.

Que sont les macros en C ?

Les macros en C sont créées avec la directive #define. Leur utilisation la plus simple est la suivante :

#define UN_NOMBRE_IMPORTANT 42

Ce qui remplace toutes les occurrences de UN_NOMBRE_IMPORTANT par la valeur littérale 42. Vous vous demandez peut-être : En quoi est-ce différent de cette définition globale ?

const int UN_NOMBRE_IMPORTANT = 42;
int main() { ... }

En pratique, pas grand-chose. Mais c’est parce que c’est très simple. Et si nous écrivions une macro fonction ?

#define MACRO_CARRE(x) x*x
int fonction_carre(int x) { return x*x; }

Quelle est la différence entre les deux ? Cela devient clair quand on considère les macros comme un simple remplacement de texte, un peu plus complexe. Elle est traitée par le préprocesseur C, qui remplace un morceau de code par un autre au moment de la compilation. Par exemple :

int x = MACRO_CARRE(2);

est développé en :

int x = 2*2;

Alors que la fonction n’est pas développée, elle se résume à une instruction d’appel. La macro présente quelques avantages. Premièrement, elle n’a pas la surcharge d’un appel de fonction. Deuxièmement, elle n’est pas typée, donc MACRO_CARRE peut être utilisée avec n’importe quel type qui supporte la multiplication.

La nature des macros qui procure ces avantages crée aussi un certain nombre de pièges dangereux. Voyons ce que MACRO_CARRE donne dans quelques cas intéressants.

MACRO_CARRE(2 + 3)

2 + 3*2 + 3

Nous voulions 25, mais nous avons obtenu 11. Pas bon. Nous pouvons facilement résoudre ce problème.

#define MACRO_CARRE (x)*(x)

Un autre exemple :

int x = 1, somme = 0;
while (x < 10) { somme += MACRO_CARRE(x++); }

est développé en :

int x = 1, somme = 0;
while (x < 10) { somme += (x++)*(x++); }

Nous voyons que x est incrémenté deux fois au lieu d’une ! En plus de cela, cela calcule x*(x+1), et non x*x.

Nous voyons donc que cette macro simple produit un résultat complètement absurde pour certaines entrées, que le compilateur considère pourtant comme valides. Bien que le fait que les macros dupliquent le code cause ces comportements dangereux, cela nous permet aussi d’être plus efficaces.

Objectifs pour la bibliothèque

Je veux que la bibliothèque soit minimaliste (header only), mais complète. Nous écrirons des fonctions pour une gamme d’opérations entre deux matrices $A$ et $B$ .

J’ai classé les fonctions nécessaires selon les critères suivants.

Type d’allocation : Devons-nous placer le résultat dans un tampon $(C := A \circ B)$ , ou l’écrire dans le premier argument $(A = A \circ B)$ ?
Type de boucle : Effectuons-nous une boucle de type produit scalaire, ou une opération élément par élément ? Devons-nous transposer le second argument $(A \circ B^{T})$ ?
Type d’opération élémentaire : Multiplie-t-on ? Additionne-t-on ? Soustrait-on ?

Notez que toutes les fonctions (sauf le produit scalaire) sont génériques par rapport à l’opération $\circ$ .

Boilerplate

Avant de commencer, nous devons définir notre type de matrice. Je vais également définir un type appelé mfloat, que nous pourrons substituer plus tard par n’importe quel type flottant.

#define DEBUG 1

typedef double mfloat;

typedef struct {
  mfloat *buf;
  int rows;
  int cols;
} matrix;

Nous allons supposer que les éléments sont stockés dans buf en ordre row-major. Maintenant, nous avons besoin de quelques fonctions basiques de type getter-setter avec une vérification optionnelle des limites, que nous utiliserons dans le reste de la bibliothèque.

static inline mfloat matrix_get(matrix m, int row, int col) {
  if (DEBUG)
    if (!(row >= 0 && col >= 0 && row < m.rows && col < m.cols)) {
      fprintf(
          stderr,
          "matrix_get: Index hors limites (%d, %d) pour une matrice "
          "de taille (%d, %d)\n",
          row, col, m.rows, m.cols);
      exit(1);
    }

  return m.buf[row * m.cols + col];
}

static inline void matrix_set(matrix m, int row, int col,
                              mfloat val) {
  if (DEBUG)
    if (!(row >= 0 && col >= 0 && row < m.rows && col < m.cols)) {
      fprintf(
          stderr,
          "matrix_set: Index hors limites (%d, %d) pour une matrice "
          "de taille (%d, %d)\n",
          row, col, m.rows, m.cols);
      exit(1);
    }

  m.buf[row * m.cols + col] = val;
}

Et nous avons besoin d’un moyen d’afficher la matrice sur la console

static inline void matrix_print(matrix m) {
  for (int i = 0; i < m.rows; i++) {
    printf("[ ");
    for (int j = 0; j < m.cols; j++) {
      // Notation scientifique, arrondie à 4 décimales
      printf("%.4e", matrix_get(m, i, j));
      printf(" ");
    }
    printf("]\n");
  }
  printf("\n");
}

Et un moyen d’allouer et de libérer des matrices sur, et depuis, le tas

matrix matrix_new(int rows, int cols) {
  double *buf = calloc(rows * cols, sizeof(double));
  if (buf == NULL) {
    printf("matrix_new: échec de calloc.");
    exit(1);
  }

  return (matrix){
      .buf = buf,
      .rows = rows,
      .cols = cols,
  };
}

## Produit scalaire

Commençons par le produit scalaire. Je vais uniquement écrire une version avec allocation de tampon, car il n'est possible d'écrire le résultat dans le premier argument que si et seulement si les deux matrices sont carrées, ce que nous ne pouvons pas supposer.

```c
static inline void matrix_dot(matrix out, const matrix m1,
                              const matrix m2) {
  if (DEBUG)
    if (m1.cols != m2.rows) {
      printf(
          "matrix dot: dimension error (%d, %d) not compat w/ "
          "(%d, %d)\n",
          m1.rows, m1.cols, m2.rows, m2.cols);
      exit(1);
    }
  for (int row = 0; row < m1.rows; row++) {
    for (int col = 0; col < m2.cols; col++) {
      double sum = 0.0;
      for (int k = 0; k < m1.cols; k++) {
        double x1 = matrix_get(m1, row, k);
        double x2 = matrix_get(m2, k, col);
        sum += x1 * x2;
      }
      matrix_set(out, row, col, sum);
    }
  }
}

## Opérations élément par élément

Chaque opération élément par élément effectue les étapes suivantes :

1. Vérifier que les deux matrices ont les mêmes dimensions
2. Exécuter l'opération sur chaque élément correspondant
3. Placer le résultat dans la matrice de sortie

On constate que la boucle sera identique pour chacune des fonctions, écrivons donc une macro pour cela.

```c
#define MAT_ELEMENTWISE_LOOP        \
  for (int i = 0; i < m1.rows; i++) \
    for (int j = 0; j < m1.cols; j++)

Et une fonction pour vérifier les dimensions, qui provoque simplement une panique si elles ne correspondent pas.

static inline void mat_bounds_check_elementwise(const matrix out,
                                                const matrix m1,
                                                const matrix m2) {
  if (DEBUG)
    if (m1.rows != m2.rows || m1.cols != m2.cols ||
        out.rows != m1.rows || out.cols != m1.cols) {
      fprintf(stderr,
              "Dimensions incompatibles pour l'opération élément par élément "
              "(%d, %d) & (%d, %d) => (%d, %d) \n",
              m1.rows, m1.cols, m2.rows, m2.cols, out.rows,
              out.cols);
      exit(1);
    }
}

Maintenant, nous voulons implémenter l’addition, la multiplication, la division et la soustraction. Comme tout le code, excepté le calcul proprement dit, est identique, nous pouvons l’abstraire dans une macro qui définit la fonction.

#define DEF_MAT_ELEMENTWISE_BUF(opname, op)           \
  static inline void matrix_##opname(                 \
      matrix out, const matrix m1, const matrix m2) { \
    mat_bounds_check_elementwise(out, m1, m2);        \
    MAT_ELEMENTWISE_LOOP {                            \
      mfloat x = matrix_get(m1, i, j);                \
      mfloat y = matrix_get(m2, i, j);                \
      matrix_set(out, i, j, op);                      \
    }                                                 \
  }

##opname insère la valeur de opname dans le nom de la fonction.

En anticipant, nous savons que nous devrons définir les fonctions d’addition, de multiplication, de division et de soustraction pour toutes les variations, alors écrivons une macro qui fait cela pour une macro de définition de fonction donnée.

#define DEF_ALL_OPS(OP_MACRO) \
  OP_MACRO(sub, (x - y));     \
  OP_MACRO(add, (x + y));     \
  OP_MACRO(div, (x / y));     \
  OP_MACRO(mul, (x * y));

Maintenant, nous pouvons réellement définir les 4 fonctions !

DEF_ALL_OPS(DEF_MAT_ELEMENTWISE_BUF)

Boum ! En une ligne, nous avons défini les fonctions matrix_add, matrix_sub, matrix_div et matrix_mul.

Essayons maintenant d’implémenter les opérations en place.

static inline void mat_bounds_check_elementwise_ip(
    matrix m1, const matrix m2) {
  if (DEBUG)
    if (m1.rows != m2.rows || m1.cols != m2.cols) {
      fprintf(stderr,
              "Dimensions incompatibles pour l'opération élément par élément en place "
              "(%d, %d) & (%d, %d) \n",
              m1.rows, m1.cols, m2.rows, m2.cols);
      exit(1);
    }
}

#define DEF_MAT_ELEMENTWISE_IP(opname, op)                 \
  static inline void matrix_ip_##opname(matrix m1,         \
                                        const matrix m2) { \
    mat_bounds_check_elementwise_ip(m1, m2);               \
    MAT_ELEMENTWISE_LOOP {                                 \
      mfloat x = matrix_get(m1, i, j);                     \
      mfloat y = matrix_get(m2, i, j);                     \
      matrix_set(m1, i, j, op);                            \
    }                                                      \
  }

Avec cette nouvelle macro, nous pouvons simplement faire

DEF_ALL_OPS(DEF_MAT_ELEMENTWISE_IP)

et matrix_ip_add, matrix_ip_mul, etc. sont définies. Nous pourrions aussi faire cela pour les opérations de transposition, mais je ne le montrerai pas ici.

Opérations unaires

Parfois, nous voulons effectuer une opération unaire sur une matrice $A$ , comme la mettre au carré $A := A^{2}$ ou la négativer $A := - A$ . Créons une fonction unaire générique pour l’opération.

#define DEF_MAT_UNARY_IP(opname, op)            \
  static inline void matrix_ip_##opname(matrix m1) { \
    MAT_ELEMENTWISE_LOOP {                           \
      mfloat x = matrix_get(m1, i, j);               \
      matrix_set(m1, i, j, op);                      \
    }                                                \
  }

DEF_MAT_UNARY_IP(square, (x * x))
DEF_MAT_UNARY_IP(negate, (-x))
DEF_MAT_UNARY_IP(sqrt, (sqrt(x)))

Et maintenant, matrix_ip_square(A), matrix_ip_negate(A), etc. sont définies.

Conclusion

Et voilà, nous avons “écrit” une bibliothèque complète d’opérations matricielles avec un effort minimal. Mais cela soulève une question : ce code est-il sûr ?

Tant que vous utilisez #undef pour tous les macros que vous avez définis, oui, il est tout aussi sûr que d’écrire toutes les fonctions vous-même. En revanche, si vous exposez des macros dans le cadre des fonctionnalités de votre bibliothèque, cela peut ne plus être sûr. Mais ce n’est pas parce qu’un code est sûr qu’il faut forcément l’utiliser. S’il y a un problème dans le code, il peut être difficile à déboguer puisque vous ne pouvez pas voir ce en quoi il est développé. Donc, dans un environnement de production, l’utilisation de macros est fortement déconseillée. Je l’ai simplement utilisée ici parce que c’était une solution de facilité pour une série éducative.

Si vous avez des questions ou des suggestions, n’hésitez pas à laisser un commentaire ou à m’envoyer un e-mail. Merci de votre lecture.

✦ Aucune IA n'a été utilisée dans la conception, la recherche, la rédaction ou l'édition de cet article.