The R Project SVN R-packages

.First.lib <- function(lib, pkg) {
    library.dynam("adapt", pkg, lib)
}


adapt <- function (ndim, lower, upper, minpts = 100, maxpts = NULL,
                   functn, eps = 0.01, xtra.lenwrk = 0:0, ...)
{
    keep.trying <- is.null(maxpts)
    minpts <- as.integer(minpts)

    if (ndim == 1) { ## fudge for 1-d functions
    warning("Using integrate() from base package for 1-d integration")
        if (keep.trying) maxpts <- minpts
    return(integrate(functn,lower,upper,subdivisions=maxpts,rel.tol=eps,...))
    }
    ## else ndim >= 2 :

    ## Check to make sure that upper and lower are reasonable lengths
    ## Both the upper and lower limits should be at least of length ndim
    if (length(lower) < ndim || length(upper) < ndim)#MM: dropped 'at least':
    stop("The lower and upper vectors need to have ndim elements\n",
         "Your parameters are:  ndim ", ndim, ", length(lower) ",
         length(lower), ", length(upper) ", length(upper))
    ff <-
    if(length(list(...)) && length(formals(functn)) > 1)
        function(x) functn(x, ...)
    else functn # .Alias

    rulcls <- 2^ndim + 2*ndim^2 + 6*ndim + 1 #-> ../src/adapt.f

    ## maxpts should be large enough.  Prefer 10*rulclc, but use 2*rulclc.
    if (keep.trying)
    maxpts <- max(minpts, 500, 2 * rulcls)
    else {
    if (minpts >= maxpts) {
        warning("maxpts must be > minpts.\n",
            "Maxpts has be increased to  minpts + 1")
        maxpts <- minpts + 1
    }
    ##
    if (maxpts < 2 * rulcls) {
        warning("You have maxpts (= ", maxpts, ") too small.\n It needs to",
            "be at least 2 times 2^ndim + 2*ndim^2 + 6*ndim+1\n",
            " and has been reset to ", 2 * rulcls, "\n")
            maxpts <- 2 * rulcls
        }
    }

    rgnstr <- 2*ndim + 3 # as in Fortran
    it <- 0
    repeat {
        it <- it+1
        ## lenwrk= length work array{allocated in C} mandated in adapt source:
    lenwrk <- as.integer(rgnstr * (1 + maxpts/rulcls)/2)
        if(verbose) cat(it,": lenwrk.orig=", lenwrk,", ")
        lenwrk <- as.integer(lenwrk + xtra.lenwrk)

    x <- .C("cadapt",
        as.integer(ndim),
        as.double(lower),
        as.double(upper),
        minpts = minpts,
        maxpts = as.integer(maxpts),
        ## now pass ff and current environment
        ff, rho = environment(),
        eps = as.double(eps),
        relerr = double(1),
        lenwrk = lenwrk,
        value = double(1),  # will contain the value of the integral
        ifail = integer(1),
                PACKAGE = "adapt")[
                c("value","relerr","minpts", "lenwrk", "ifail", "eps")]

        if(verbose)
            cat("  needed (lenwrk equivalent of)", x$eps+ 2*(ndim + 3)+ 2,"\n")
    if (x$ifail == 1 && keep.trying)
        maxpts <- maxpts*2
    else
        break
    }
    if(x$ifail)
    warning(x$warn <-
        c("Ifail=1, maxpts was too small. Check the returned relerr!",
          paste("Ifail=2, lenwrk=",lenwrk," was too small, needed",
                        ## sbrgns + rgnstr + 2 (<- ../src/adapt.f )  =
                        x$eps + 2*(ndim + 3) + 2,
                        ".\n ***> fix adapt(); check the returned relerr!"),
          "Ifail=3: ndim > 20 -- rewrite the fortran code ;-) !",
          "Ifail=4, minpts > maxpts; should not happen!",
          "Ifail=5, internal non-convergence; should not happen!",
          )[x$ifail])

    x$eps <- NULL
    x$iter <- it
    class(x) <- "integration"
    x
}

print.integration <- function(x, ...) {
    print(noquote(sapply(x, format, ...)),...)
    invisible(x)
}