[R] Calling Rdqags doesn't produce correct result.

Paul August paulaugust2003 at yahoo.com
Thu Nov 18 04:55:59 CET 2004


Does anyone has a clue what went wrong in the
following attempt?

I am trying to call the R built-in function Rdqags()
from my C
program for numerical integration. Following are the C
program
and the corresponding R program:

C program
---------
    void test(double *a,
              double *b,
              double *epsabs,
              double *epsrel,
              double *result,
              double *abserr,
              int *neval,
              int *ier,
              int *limit,
              int *lenw,
              int *last,
              int *iwork,
              double *work,
              double *exx)
    {
        void *ex;
    
        ex = exx;
        Rdqags(tmpfun, ex, a, b, epsabs, epsrel,
result, abserr, neval, ier,
                limit, lenw, last, iwork, work);
    }
    
    // User supplied function
    void tmpfun(double *x, int n, void *ex)
    {
        int i;
        double *tmp;
    
        tmp = (double *)ex;
    
    //    for(i=1;i<=n;i++) {x[i] = pow(x[i], *tmp);}
        for(i=1;i<=n;i++) {x[i] = 2.0*x[i];}
        return;
    }



R program
---------
    test <- function(a, b, epsabs =
.Machine$double.eps^0.25, epsrel = epsabs,
                limit = 100, lenw = 400, ex = -0.5)
    {
        val <- .C("test",
            as.double(a),
            as.double(b),
            as.double(epsabs),
            as.double(epsrel),
            result = double(1),
            abserr = as.double(-1),
            neval = as.integer(-1),
            ier = as.integer(-1),
            as.integer(limit),
            as.integer(lenw),
            last = as.integer(-1),
            integer(limit+1),
            double(lenw+1),
            as.double(ex),
            NAOK=T,
            specialsok=T
        )
        val
    }


I followed the instructions in Writing R Extensions ->
The R API
-> Integration, and the program was complied and
loaded into R
successfully. However, calling R function test(0,
20)$result
gives me 385.0554, not 400! The returned ier value of
test(0, 20)
is 0! But if I call integrate(function(x) 2 * x, lower
= 0, upper
= 20), I got the correct result of 400. Apparently,
there must be
something wrong in my programs. But I cannot see where
is the
mistake. Do I miss anything obvious here?

The programs above were tested in R-1.9.0 for
MS-Windows and the
compiler is MSC++ V6. I have used similar calls to the
built-in
random generators and never had any problems.

Paul.




More information about the R-help mailing list