Appel d'un programme interne à R dans un programme C

Auteur.e.s
Annie Bouvier, Inra
Résumé

Un certain nombre des programmes internes de R sont rendus facilement utilisables dans un programme C: ce sont les programmes de l’API. Ils effectuent diverses tâches courantes, dont l’impression, l’allocation mémoire, la gestion des erreurs, la génération de nombres aléatoires, l’optimisation, l’intégration.

Notations: RHOME: chemin du répertoire d’installation de R sur le site de l’utilisateur (pour afficher celui-ci, tapez: R RHOME). Le terme “programme” est ici utilisé pour désigner un programme ou une “subroutine” C ou fortran. On réserve le terme “fonction” pour ce qui est du code R.

Appel d’un programme de l’API R pour C

Un certain nombre des programmes internes de R sont rendus facilement utilisables dans un programme C: ce sont les programmes de l’API. Ils effectuent diverses tâches courantes, dont l’impression, l’allocation mémoire, la gestion des erreurs, la génération de nombres aléatoires, l’optimisation, l’intégration. L’API contient aussi des constantes mathématiques (PI, Inf, etc …).

Lorsque le programme C est appelé depuis R, il est vivement conseillé d’utiliser les programmes de l’API, de préférence à ceux du compilateur ou aux siens propres, dès lors qu’ils ont une répercussion sur le système pour éviter toute interférence avec R (par exemple l’impression ou l’allocation mémoire).

Exemple: appel d’un programme de génération de nombres aléatoires

Le programme C, contenu dans le fichier ma_norm.c appelle le programme de l’API norm_rand générateur de nombres aléatoires suivant une loi normale:

// inclusion du fichier de R qui déclare les programmes
#include <R.h> 
void ma_norm(int *N, double *res) {
  int i;
  GetRNGstate(); // lecture de la graine du générateur
  for (i=0; i< *N; i++) {
    res[i] = norm_rand();
//  Rprintf est un programme de l'API: il s'utilise comme printf
    Rprintf("res[%d] = %f\n", i, res[i]); 
  } // fin i
  PutRNGstate(); // écriture de la graine du générateur
}

la compilation/édition de liens est réalisée avec la commande suivante et crée la librairie ma_norm.so:

R CMD SHLIB ma_norm.c

Utilisation.

dyn.load("ma_norm.so")  # chargement de la librairie
n <- 5
res <- vector(mode = "numeric", length = n)  # allocation des sorties
out <- .C("ma_norm", as.integer(n), res = as.double(res))$res
print(out)
res[0] = -0.947863
res[1] = 0.158358
res[2] = -0.316687
res[3] = -0.527797
res[4] = -0.145287
[1] -0.9478626  0.1583578 -0.3166872 -0.5277972 -0.1452874

Appel d’un programme de l’API lorsque le programme C n’est pas appelé depuis R

Si l’installation de R le permet, c’est-à-dire si cela a été prêvu à l’étape de configuration, l’utilisation de la librairie Rmath.h de R est possible depuis un programme C autonome. L’exemple précédent s’écrirait ainsi dans un programme C autonome:

// inclusion du fichier de R qui déclare les programmes
#define MATHLIB_STANDALONE
#include  <Rmath.h>
#include  <stdio.h>
void ma_norm(int N, double *res) {
  int i;
  for (i=0; i< N; i++) {
    res[i] = norm_rand();

    printf("res[%d] = %f\n", i, res[i]); 
  } // fin i

}
int main() {
  int N=5;
  double res[N];
  ma_norm(N, res);
  return(0);
}

La commande de compilation doit spécifier l’emplacement de la librairie Rmath qui contient les programmes R. Cet emplacement varie selon les installations.

Par exemple:

cc ma_norm.c -L/usr/local/lib/R/lib -lRmath -lm

Appel d’un programme d’une librairie de R

R utilise pour ses propres besoins les librairies d’algèbre linéaire BLAS, LAPACK et LINPACK/EISPACK. Les fichiers de RHOME, R_ext/BLAS.h, R_ext/Lapack.h et R_ext/Linpack.h contiennent les déclarations des programmes de ces librairies qui sont à la disposition de l’utilisateur. Ce sont toutes des points d’entrée Fortran.

Exemple: appel d’un programme Fortran de LAPACK dans un programme C

Nous utilisons le programme DGESDD de LAPACK (singular-value decomposition of a rectangular matrix).

Le programme C, contenu dans le fichier monsvd.c.

#include <R.h>
#include <Rinternals.h>
#include <Rdefines.h>
#include <R_ext/Lapack.h>
#include <R_ext/Utils.h>
SEXP monsvd(SEXP x, SEXP res) {
  /* Déclarations */
 int *xdims, n, p, lwork;
 double u=1, v=1;
 int  ldu=1, ldvt=1;
 int info =0;
 double *work, *xvals, tmp;
 char jobu='N';

  /* Accès aux dimensions de x */
 xdims = INTEGER(coerceVector(getAttrib(x, R_DimSymbol), INTSXP));
 n = xdims[0]; p = xdims[1];

/* Allocation d'un vecteur de travail double */
 xvals = (double *) R_alloc(n * p, sizeof(double));
    /* work on a copy of x */
    Memcpy(xvals, REAL(x), n * p);
  {
    /* Allocation d'un vecteur de travail entier */
    int *iwork= (int *) R_alloc(8*(n<p ? n : p), sizeof(int));

    /* ask for optimal size of work array */
    lwork = -1;
    F77_CALL(dgesdd)(&jobu,
             &n, &p, xvals, &n, REAL(res),
             &u, &ldu,
             &v, &ldvt,
             &tmp, &lwork, iwork, &info);
    if (info != 0)
        error("error code %d from Lapack routine '%s'", info, "dgesdd");
    lwork = (int) tmp;
    work = (double *) R_alloc(lwork, sizeof(double));
    F77_CALL(dgesdd)(&jobu,
             &n, &p, xvals, &n, REAL(res),
             &u, &ldu,
             &v, &ldvt,
             work, &lwork, iwork, &info);
    if (info != 0)
        error("error code %d from Lapack routine '%s'", info, "dgesdd");
    }
  return(res);
}

Remarque: Les arguments monsvd sont déclarés de type SEXP, ce qui rend disponible certaines de leurs caractéristiques, comme ici les dimensions de x.

Compilation/édition de liens.

Nous devons spécifier où se trouve la librairie LAPACK lors de l’édition de liens, ce qui se fait en créant un fichier nommé Makevars dans le directory des sources, fichier que la commande R CMD SHLIB prend automatiquement en compte (attention: il faut le nommer obligatoirement Makevars).

Ce fichier sert à définir des variables de compilation et d’édition de liens. La variable PKG_LIBS définit les librairies. Notre fichier Makevars contient une seule ligne:

PKG_LIBS=$(LAPACK_LIBS)

$(LAPACK_LIBS) est une variable reconnue par R CMD SHLIB: elle est égale à la localisation de la librairie LAPACK dans l’installation courante de R.

Nous compilons par la commande usuelle:

R CMD SHLIB  monsvd.c

Utilisation.

dyn.load("monsvd.so")  # Chargement de la bibliothèque
n <- 9
i <- 1:n
X <- 1/outer(i - 1, i, "+")[, 1:6]
u <- rep(-9, ncol(X))
# Appel du programme
print(.Call("monsvd", X, u))
[1] 1.668433e+00 2.773727e-01 2.223722e-02 1.084693e-03 3.243788e-05
[6] 5.234864e-07

Appel d’un programme interne à un paquetage de R

L’exemple suivant illustre l’appel à un programme utilisé dans un paquetage de R, en l’occurrence le programme inplg du paquetage sgeostat, qui permet de savoir si un ou plusieurs points d’un espace à 2 dimensions sont à l’intérieur d’un polygone. Il est utilisé dans la fonction in.polygon(), que nous affichons pour connaitre la liste des arguments:

in.polygon
    function (x0, y0, x, y) 
{
    n <- length(x)
    n0 <- length(x0)
    if (length(y0) != n0) 
        stop("length of x0 and y0 differ!")
    ret <- .Fortran("inplg", as.double(c(x, 0)), as.double(c(y, 
        0)), as.integer(n), as.double(x0), as.double(y0), as.integer(n0), 
        inhull = integer(n0))
    as.logical(ret$inhull)
}

Le programme Fortran, contenu dans le fichier appelf.f. Le programme inplg étant en Fortran, nous écrivons le programme qui l’appelle dans le même langage:

          SUBROUTINE APPELF(XP,YP,NP,X0,Y0,N0, INHULL )
*     .. Array Arguments ..
      DOUBLE PRECISION                 XP(NP), YP(NP)
      DOUBLE PRECISION                 X0(N0),Y0(N0)
      INTEGER                          INHULL(N0)
*     .. Scalar Arguments ..
      INTEGER                          NP,N0
*     .. External Subroutines ..
      EXTERNAL INPLG
* -----------------------------------------------
*     .. Local Scalars ..
      INTEGER   I

      CALL INPLG(XP,YP, NP, X0, Y0, N0, INHULL)
         END

Compilation/édition de liens.

Comme dans l’exemple précédent, un fichier Makevars indique où se trouve la librairie qui contient le programme inplg:

    PKG_LIBS= RHOME/library/sgeostat/libs/sgeostat.so

Nous compilons par la commande usuelle:

R CMD SHLIB appelf.f

Utilisation.

dyn.load("appelf.so")  # Chargement de la bibliothèque
# Création du polygone:
x <- c(0, 1, 0, -1, 0)
y <- c(-1, 0, 1, 0, 0)
n <- length(x)
# Les points à tester:
x1 <- c(0, 1, 0.5, 0.5)
y1 <- c(0, 1, 1, 0)
n1 <- length(x1)
# Allocation des sorties:
inhull <- vector(mode = "numeric", length = n1)
# Appel du programme Fortran
res <- .Fortran("appelf", as.double(x), as.double(y), as.integer(n), as.double(x1), 
    as.double(y1), as.integer(n1), inhull = as.integer(inhull))$inhull
print(res)
[1] 1 0 0 1

Références

La liste des programmes de l’API figure dans la documentation de R: Writing Extensions au chapitre “R API: entry points for C code”. Contrairement à ce que le titre pourrait suggérer, des points d’entrée Fortran y sont aussi présentés.

Le source des librairies d’algèbre linéaire n’est pas disponible dans une installation standard de R. On les trouvera en déchargeant le code source de R depuis le site du CRAN. Ils sont dans le répertoire src/modules.

Versions des outils utilisés
R version 3.4.2 (2017-09-28)
Thèmes de la fiche
Thèmes