2013-08-28 9 views
5

Napisałem kawałek kodu C, który deklaruje kwadratową macierz rozmiaru 4x4. Następnie pobiera próbki z funkcji próbkowania o nazwie rgig w opakowaniu GeneralizedHyperbolic w R. Odwraca ona matrycę przy użyciu biblioteki gsl od gnu i wypluwa wynik. Jest to ćwiczenie w nazywaniu R od Cz wbudowanym R w C

#include <stdio.h> 
#include <stdio.h> 
#include <math.h> 
#include <stdlib.h> 
#include <stddef.h> 

// for gsl 
#include <gsl/gsl_machine.h> 
#include <gsl/gsl_rng.h> 
#include <gsl/gsl_randist.h> 
#include <gsl/gsl_cdf.h> 
#include <gsl/gsl_cblas.h> 
#include <gsl/gsl_sf_gamma.h> 
#include <gsl/gsl_vector.h> 
#include <gsl/gsl_matrix.h> 
#include <gsl/gsl_blas.h> 
#include <gsl/gsl_linalg.h> 

// for R embedding in C 
#include <Rinternals.h> 
#include <Rdefines.h> 
#include <Rembedded.h> 
#include <R_ext/Parse.h> 

void gsl_square_matrix_inverse (gsl_matrix *, gsl_matrix *, int); 
SEXP get_rInvGauss(void); 

int main(void) 
{ 
    // Define the dimension n of the matrix 
    // and the signum s (for LU decomposition) 
    int s, i, j, n = 4; 
    // Define all the used matrices 
    gsl_matrix * m = gsl_matrix_alloc (n, n); 
    gsl_matrix * inverse = gsl_matrix_alloc (n, n); 

    // R embedding in C 
    char *localArgs[] = {"R", "--no-save","--silent"}; 
    SEXP rInvGauss; 

    // init R embedding 
    Rf_initEmbeddedR(3, localArgs); 

    printf("\n Printing matrix m before set. size %d by %d... \n", n, n); 
    for(i=0; i<n; i++){ 
     printf("\n"); 
     for(j=0; j<n; j++){ 
      printf(" %f ", gsl_matrix_get(m, i, j)); 
     } 
    } 

    // set diagonal elements of matrix m from Inverse Gaussian Random samples 
    for(i=0; i<n; i++){ 
     rInvGauss = get_rInvGauss(); 
     gsl_matrix_set(m, i, i, *REAL(rInvGauss)); 
    } 

    Rf_endEmbeddedR(0); // end the R embedding in C 

    printf("\n Printing matrix m ..... \n"); 
    for(i=0; i<n; i++){ 
     printf("\n"); 
     for(j=0; j<n; j++){ 
      printf(" %f ", gsl_matrix_get(m, i, j)); 
     } 
    } 

    // inverse of matrix m 
    gsl_square_matrix_inverse (m, inverse, n); 

    printf("\n Printing inverse of matrix m ..... \n"); 
    for(i=0; i<n; i++){ 
     printf("\n"); 
     for(j=0; j<n; j++){ 
      printf(" %f", gsl_matrix_get(inverse, i, j)); 
     } 
    } 

    return 0; 
} 


SEXP get_rInvGauss(void) { 
    SEXP e, s, t, tmp, result; 
    int errorOccurred, n=1; 
    double chi=5, psi=4, lambda=0.5; 

    // create and evaluate 'require(GeneralizedHyperbolic)' 
    PROTECT(e = lang2(install("require"), mkString("GeneralizedHyperbolic"))); 
    R_tryEval(e, R_GlobalEnv, &errorOccurred); 
    if (errorOccurred) { 
     // handle error 
     printf("\n Error loading library GeneralizedHyperbolic:"); 
    } 
    UNPROTECT(1); 


    // Create the R expressions using a paired list 
    // rgig(n = 1, chi = 5, psi = 4, lambda = 0.5) with the R API. 
    PROTECT(t = s = allocVector(LANGSXP, 5)); 
    // could also be done by: PROTECT(t = s = allocList(5)); SET_TYPEOF(s, LANGSXP); 

    tmp = findFun(install("rgig"), R_GlobalEnv); 
    if(tmp == R_NilValue) { 
     printf("No definition for function rgig.\n"); 
     UNPROTECT(1); 
     exit(1); 
     } 
    SETCAR(t, tmp); t = CDR(t); 
    SETCAR(t, ScalarInteger(n)); SET_TAG(t, install("n")); t= CDR(t); 
    SETCAR(t, ScalarReal(chi)); SET_TAG(t, install("chi")); t= CDR(t); 
    SETCAR(t, ScalarReal(psi)); SET_TAG(t, install("psi")); t= CDR(t); 
    SETCAR(t, ScalarReal(lambda)); SET_TAG(t, install("lambda")); t= CDR(t); 
    PROTECT(result = R_tryEval(VECTOR_ELT(s, 0), R_GlobalEnv, NULL)); 
    UNPROTECT(2); 

    return(result); 

} 

void gsl_square_matrix_inverse (gsl_matrix *m, gsl_matrix *inverse, int n){ 

    int s, i, j; 
    gsl_permutation * perm = gsl_permutation_alloc (n); 
    // Make LU decomposition of matrix m 
    gsl_linalg_LU_decomp (m, perm, &s); 
    // Invert the matrix m 
    gsl_linalg_LU_invert (m, perm, inverse); 

} 

skompilowany kod za pomocą:

R CMD SHLIB -lgsl -lgslcblas embedR_matinv.c 

z wyjściem:

gcc -arch x86_64 -std=gnu99 -I/Library/Frameworks/R.framework/Resources/include -I/Library/Frameworks/R.framework/Resources/include/x86_64 -DNDEBUG -I/usr/local/include -fPIC -g -O2 -c embedR_matinv.c -o embedR_matinv.o 
gcc -arch x86_64 -std=gnu99 -dynamiclib -Wl,-headerpad_max_install_names -undefined dynamic_lookup -single_module -multiply_defined suppress -L/usr/local/lib -o embedR_matinv.so embedR_matinv.o -lgsl -lgslcblas -F/Library/Frameworks/R.framework/.. -framework R -Wl,-framework -Wl,CoreFoundation 

Kiedy złożyć:

R CMD embedR_matinv 

its daje błąd:

/Library/Frameworks/R.framework/Resources/bin/Rcmd: line 62: exec: embedR_matinv: not found 

Co robię źle?

Ja również zmienił main() do test() i wykonał obiekt udostępniony jako

R CMD SHLIB -lgsl -lgslcblas embedR_matinv.c -o embedR_matinv 

z wyjściem:

gcc -arch x86_64 -std=gnu99 -I/Library/Frameworks/R.framework/Resources/include -I/Library/Frameworks/R.framework/Resources/include/x86_64 -DNDEBUG -I/usr/local/include -fPIC -g -O2 -c embedR_matinv.c -o embedR_matinv.o 

Jeśli robię dyn.load("embedR_matinv.so") w R Studio, kod działa bez wypowiedzenia tzn. wisi!

Wszelkie sugestie na temat tego, co jest nie tak w kodzie?

+0

Co jest generowany podczas kompilacji 'embedR_matinv.c'? Czy jesteś pewien, że musisz wydać 'R CMD embedR_matinv', a nie' R CMD embedR_matinv.o' lub 'R CMD embedR_matinv.so'? Może powinieneś podać pełną ścieżkę do wygenerowanego pliku? – mangusta

+0

Gdzie w kodzie jest zawieszony? Pachnie, jakby utknął w tłumaczu R we wbudowanym kodzie R. –

Odpowiedz

0

Czy używasz systemu Mac OSX?

spróbuj tego:

sudo cp /Library/Frameworks/R.framework/Resources/library/Rserve/libs/x86_64/Rserve-bin.so\ /Library/Frameworks/R.framework/Resources/bin/Rserve