1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
|
## debugging example
success <- c(13,12,11,14,14,11,13,11,12)
failure <- c(0,0,0,0,0,0,0,2,2)
resp <- cbind(success, failure)
predictor <- c(0, 5^(0:7))
### try would suppress a traceback
# glm(resp ~ 0+predictor, family = binomial(link="log"))
# traceback()
## R code to run the .Call/.External examples in
## `Writing R extensions'
dyn.load(paste("R-exts", .Platform$dynlib.ext, sep=""))
## ----- outer products example -----
out <- function(x, y)
{
storage.mode(x) <- storage.mode(y) <- "double"
.Call("out", x, y)
}
out(1:3, 2:4)
x <- 1:3; names(x) <- letters[x]
out(x, 2:4)
## ----- convolution example -----
conv <- function(a, b) .Call("convolve2", a, b)
u <- rep(1, 5)
conv(u, u)
convE <- function(a, b) .External("convolveE", a, b)
convE(u, u)
## ----- Lists examples -----
showArgs <- function(...) invisible(.External("showArgs", ...))
showArgs(u = u, x = x, let = letters)
showArgs1 <- function(...) invisible(.Call("showArgs1", list(...)))
showArgs1(u = u, x = x, let = letters)
a <- list(a = 1:5, b = rnorm(10), test = runif(100))
.Call("lapply", a, quote(sum(x)), new.env())
.Call("lapply2", a, sum, new.env())
## ----- zero-finding -----
zero <- function(f, guesses, tol = 1e-7) {
f.check <- function(x) {
x <- f(x)
if(!is.numeric(x)) stop("Need a numeric result")
as.double(x)
}
.Call("zero", body(f.check), as.double(guesses), as.double(tol),
new.env())
}
cube1 <- function(x) (x^2 + 1) * (x - 1.5)
zero(cube1, c(0, 5))
## ----- numerical derivatives -----
numeric.deriv <- function(expr, theta, rho = sys.frame(sys.parent()))
{
eps <- sqrt(.Machine$double.eps)
ans <- eval(substitute(expr), rho)
grad <- matrix(, length(ans), length(theta),
dimnames = list(NULL, theta))
for (i in seq_along(theta)) {
old <- get(theta[i], envir=rho)
delta <- eps * min(1, abs(old))
assign(theta[i], old+delta, envir=rho)
ans1 <- eval(substitute(expr), rho)
assign(theta[i], old, envir=rho)
grad[, i] <- (ans1 - ans)/delta
}
attr(ans, "gradient") <- grad
ans
}
omega <- 1:5; x <- 1; y <- 2
numeric.deriv(sin(omega*x*y), c("x", "y"))
.External("numeric_deriv", quote(sin(omega*x*y)),
c("x", "y"), .GlobalEnv)
|
