Appeler du C, du C++ ou du Fortran depuis R. Notions avancées

Vecteurs, matrices ou tableaux en argument

void remp_vect(double *vect)
{
/* Remplir un vecteur de longueur >=5 par des entiers consécutifs */
  int i;
  for (i=0; i<5; i++)
    vect[i]=i;
}

# Auparavant sous Unix, on a créé la librairie par: R CMD SHLIB
# remp_vect.c

# Chargement de la librairie
dyn.load("remp_vect.so")
vect <- rep(0, 5)  #allocation du vecteur résultat
out <- .C("remp_vect", retour = as.double(vect))
out$retour  # la valeur de retour de la fonction 

void remp_vect(double *vect, int *taille)
{
/* Remplir un vecteur de longueur >=taille d'entiers consécutifs */
  int i;
  for (i=0;i< *taille;i++)
    vect[i]=i;
}

# Auparavant sous Unix, ON A RECRÉE LA LIBRAIRIE PAR: R CMD SHLIB
# remp_vect.c

# Chargement de la librairie
dyn.load("remp_vect.so")
vect <- rep(0, 5)  #allocation du vecteur résultat
out <- .C("remp_vect", retour = as.double(vect), as.integer(length(vect)))
out$retour  # la valeur de retour de la fonction 

[1] 0 1 2 3 4

void somme_mat(double *mat, int *lignes, int *colonnes)
{
  int l,c,nl,nc;

  nl=*lignes;
  nc=*colonnes - 1;

  for (l=0;l<nl;l++) {
    /* Initialisation a zero du dernier element de la ligne 'l',
       c'est-à-dire:    mat[l, nc]=0 */
    mat[nl*nc+l]=0;
    for (c=0;c<nc;c++) 
      /* Exécution de:   mat[l, nc]+=mat[l, c] */
      mat[nl*nc+l]+=mat[nl*c+l];
  }
}

# Auparavant sous Unix, création de la librairie: R CMD SHLIB
# somme_mat.c
dyn.load("somme_mat.so")
mat <- matrix(data = as.double(c(1, 4, 7, 2, 5, 8, 3, 6, 9, 0, 0, 0)), 
    nrow = 3, ncol = 4)
out <- .C("somme_mat", retour = mat, as.integer(nrow(mat)), as.integer(ncol(mat)))
out$retour

     [,1] [,2] [,3] [,4]
[1,]    1    2    3    6
[2,]    4    5    6   15
[3,]    7    8    9   24

Définir les options de compilation via un Makevars

Par défaut, la compilation par R CMD COMPILE ou R CMD SHLIB est effectuée en utilisant les compilateurs et les options qui ont servi à l’installation de R sur la machine courante, options que l’on peut voir dans le fichier /etc/Makeconf du répertoire d’installation de R.

Pour définir ses propres options, créer un fichier Makevars dans le répertoire des sources.

Exemple (fichier Makevars):

`PKG_CPPFLAGS=-Wall

Dans cet exemple, l’option de compilation “warning bavard” est requise. Lorsque la commande R CMD SHLIB est exécutée, cette option remplace ou se rajoute aux options standards. Voir un autre exemple d’utilisation de Makevars au paragraphe Utilisation d’une librairie externe: exemple de la librairie `Nag.

Précautions à prendre dans ses programmes (imprimer, gérer les erreurs, allouer de la mémoire, autoriser les interruptions, générer des nombres aléatoires)

Pour ne pas interférer avec la gestion de la session opérée par R, il est fortement conseillé d’utiliser les programmes de l’API R pour C, lorsqu’il en existe pour la tâche à réaliser, c’est-à-dire pour imprimer, gérer les erreurs, allouer de la mémoire, autoriser les interruptions, générer des nombres aléatoires. Il faut respecter aussi certaines conventions lorsqu’on appelle du Fortran depuis le C ou l’inverse.

Exemple: Le programme C (fichier prec.c)

/* Inclusion des fichiers R qui déclarent les types et fonctions de
l'API */
#include <R.h>
#include <R_ext/Utils.h> // pour permettre d'interrompre
void prec(int *N, double *res) {
  int i;
  /* R_alloc: allocation pour la durée du programme C */ 
  double *coeffs = ( double *) R_alloc((*N), sizeof(double));
  if (!coeffs)
    error("coeffs nul");  /* error: programme d'erreur de R */
  Rprintf("Début calcul\n");/* Rprintf: programme d'écriture de R */
  GetRNGstate(); /* GetRNGstate: lire la graine du générateur de R */
 for ( i=0; i< *N; i++) {
    R_CheckUserInterrupt(); /* R_CheckUserInterrupt: permettre les interruptions */
    coeffs[i] = norm_rand(); /* norm_rand: appel du générateur de
                              nombres  aléatoires de R */
    } // fin i
  PutRNGstate(); /* PutRNGstate: écrire la graine du générateur */
  /* Effectuer des calculs
.....
  */
  /* Sommer les résultats */
  *res=0.0;
  for ( i=0; i< *N; i++)
    *res += coeffs[i];
  Rprintf("Fin calcul, res=%g\n", *res);
}

Remarque sur R_CheckUserInterrupt: Lorsque l’utilisateur tape sur Ctrl-C pour interrompre un programme interfacé à R, il ne se passe rien: l’interruption n’est pas effectuée. Pour qu’elle le soit, il faut que le programmeur l’autorise explicitement en insérant un appel à la fonction R_CheckUserInterrupt() à l’intérieur des boucles qu’il veut rendre interruptibles.

Le programme R : (fichier prec.R)

# Auparavant sous Unix, création de la librairie: R CMD SHLIB prec.c
dyn.load("prec.so")
out <- .C("prec", as.integer(5), retour = as.double(0))

Début calcul
Fin calcul, res=-1.11854

out$retour

[1] -1.118536

Voir aussi la fiche Appel d’un programme interne à R dans un programme C et la documentation R: 6. The R API: entry points for C code dans l’aide en ligne de R.

Utilisation d’une librairie externe: exemple de la librairie Nag

La bibliothèque `NAG est une bibliothèque de sous-programmes FORTRAN, d’analyse numérique et statistique. Nous décrivons sur un exemple, les étapes nécessaires à son utilisation depuis R.

1. Créer, dans le répertoire de travail, un fichier nommé Makevars. Il définit les options nécessaires à l’éditeur de liens, c’est-à-dire ici, le chemin de la librairie NAG:

PKG_LIBS=-L/usr/local/library  -lnag_nag

2. Créer le programme qui appelle NAG. Ici, le programme Fortran BFT appelle la subroutine NAG X02AJF: (fichier `bft.f)

C     A MACHINE-DEPENDENT CONSTANT IS RETURNED HERE. EPSMCH IS THE
C     SMALLEST POSITIVE REAL NUMBER SUCH THAT 1.0 + EPSMCH .GT. 1.0
C
      SUBROUTINE BFT(EPSMCH)
       DOUBLE PRECISION EPSMCH
      EPSMCH = X02AJF()
      write(6,1000) 'EPSMCH=', EPSMCH
 1000  FORMAT(A,  D15.8)
      RETURN
      END

3. Créer la librairie partagée:

R CMD SHLIB bft.f

4. Dans la session R, charger la librairie partagée:

dyn.load("bft.so")

et utiliser (fichier `bft.R`):

.Fortran("bft", as.double(0))

EPSMCH= 0.14012985D-44
[[1]]
[1] 1.401298e-45

Appel d’un programme C++

Pour appeler des programmes C++ depuis R, il faut les enrober dans une déclaration extern “C”. La compilation se fait par la commande usuelle R CMD SHLIB, qui reconnait les fichiers C++ d’après le suffixe de leurs noms.

Pour plus d’information, voir Interfacing C++ code dans l’aide en ligne de R, et le package Rcpp.

Deboguer un code C ou Fortran

 R -d gdb # Invoquer R avec les options qui lancent le mode debug
(gdb) run
for(i in 1:1e7) x <- rnorm(100)
Interrompre
(gdb) where
#0  Rf_qnorm5 (p=0.67359275277125485, mu=, sigma=1, 
lower_tail=1, log_p=0) at qnorm.c:148
etc ...

(gdb) quit

`sprof` est une commande Unix permettant de "profiler" des librairies partagées, c'est-à-dire de compter les nombres d'appels aux fonctions lors d'une exécution, de générer le graphe des appels, etc... Elle peut être utilisée pour "profiler" des librairies créées par `R CMD SHLIB`.

export LD_PROFILE_OUTPUT=/

dyn.load("somme_mat.so")
mat <- matrix(data = as.double(c(1, 4, 7, 2, 5, 8, 3, 6, 9, 0, 0, 0)), 
    nrow = 3, ncol = 4)
out <- .C("somme_mat", retour = mat, as.integer(nrow(mat)), as.integer(ncol(mat)))

sprof  -p /home/paul/projets/somme_mat.so
               /home/paul/projets/somme_mat.so.profile
Flat profile:

Each sample counts as 0,01 seconds.
 %   cumulative   self              self     total
time   seconds   seconds    calls  us/call  us/call  name

Passer une liste en argument d’un programme C; utilisation de Call

Ce paragraphe explique comment accéder à une structure complexe de R — une liste — dans un programme C. Il illustre aussi l’utilisation de la commande Call, qui, à l’inverse de .C(), ne duplique pas les arguments.

Le programme C (fichier scr.c) inclut les fichiers de R qui déclarent les types et les fonctions nécessaires (instructions #include). Nous avons ensuite écrit la fonction getListElement qui renvoie un composant de liste, sachant son nom. La fonction calcSCR l’utilise pour extraire le composant de nom lesMod de la liste modListe passée en argument. Ce composant est lui-même une liste, qui comprend au moins 3 composants. On en extrait le 3ième, qui est un vecteur de 2 nombres entiers. Ces nombres sont modifiés et la liste est retournée.

#include <R.h>
#include <Rinternals.h>
#include <Rdefines.h>
/* ++++++++++++++++++++++++++++++++++++++++++++++++
FUNCTION
 Access to a component by name from a R-list
INPUT
 list: a R-list
 str: name of a component of list
RETURN VALUE
 the component "str" of list
 ++++++++++++++++++++++++++++++++++++++++++++++++ */
SEXP getListElement(SEXP list, const char *str)
{
/* Initialiser le retour */
  SEXP elmt = R_NilValue;
 /* Accéder aux noms des composants */
  SEXP names=getAttrib(list, R_NamesSymbol);
  int i;
/* Recherche du composant de nom voulu */
  for (i=0; i< length(list); i++)
    if (strcmp(CHAR(STRING_ELT(names, i)), str) ==0) {
      elmt=VECTOR_ELT(list, i);
      break;
    }
  return(elmt);
} /* Fin de getListElement */

SEXP calcSCR( SEXP modListe) {
/* Structures de travail */
  SEXP   listlesMod ;
  int *matlesMod;
/* Acces au composant "lesMod" de la liste en argu */
  SEXP glesMod = getListElement(modListe,"lesMod");
/* Prendre le 3ieme composant, i.e celui d'indice 2
 (l'indiçage commence à 0) */
  listlesMod = VECTOR_ELT(glesMod, 2);
/* C'est un vecteur d'entiers; on accède à son pointeur */
   matlesMod = INTEGER_POINTER(listlesMod);
/* Modification de son 1ier élément */
   matlesMod[0] = matlesMod[0] +1;
/* Retourner la liste modifiée */
   return(modListe);
}

Remarques:

• Les arguments de type list sont de type SEXP, type interne à R déclaré dans les fichiers include.

• L’accès à un composant de liste d’après sa position, se fait en utilisant la fonction VECTOR_ELT.

• L’accès à des valeurs entières se fait en utilisant la fonction INTEGER_POINTER, à des valeurs réelles en utilisant la fonction NUMERIC_POINTER et à des chaînes de caractères en utilisant la fonction STRING_ELT.

• Les arguments déclarés SEXP ne sont pas dupliqués lors de l’appel au C: ils sont passés par adresse. La modification de l’argument formel `modListe implique celle de l’argument réel.

L’appel dans R: (fichier scr.R)

# Auparavant sous Unix, création de la librairie: R CMD SHLIB scr.c
dyn.load("scr.so")
lesMod <- list(a = 1, b = 2, c = as.integer(c(1, 2)))
z <- list(lesMod = lesMod)
z.out <- .Call("calcSCR", z)  # Appel par .Call
print(z.out)

$lesMod
$lesMod$a
[1] 1

$lesMod$b
[1] 2

$lesMod$c
[1] 2 2

print(z)  # z est aussi modifié

$lesMod
$lesMod$a
[1] 1

$lesMod$b
[1] 2

$lesMod$c
[1] 2 2

Références

Références

System and foreign language interfaces, dans l’aide en ligne de R, rubrique “Writing R extensions”.

Voir aussi

La fiche Appel d’un programme interne à R dans un programme C