Skip to content
Prev 277967 / 398506 Next

Function gives numeric(0) for every input

Hi Alex,

Michael is correct that vectorizing this is the way to go.  Since this
is for class and vectorizing requires a bit of thought, I am not going
to demonstrate it with your function; however, here are some very
minor changes you can make that give fairly substantial performance
increases.  It can be hard when you are learning create super
efficient code, but as a general rule of thumb (particularly when
loops are involved), byte compiling cannot hurt and may help quite a
bit, and it does not really require any extra thought or effort on
your part.

Cheers,

Josh (examples and benchmarks below)


## what you had just changing phat <- 0L
## and indenting the code and using <- instead of =
buffonslow <- function(n) {
  x <- NULL
  theta <- NULL
  a <- NULL
  phat <- 0L
  i <- 1L
  while (i <= n) {
    x[i] <- runif(1, 0, 1/2)
    theta[i] <- runif(1, 0, pi/2)
    a[i] <- 1/2 * (sin(theta[i])) - x[i]
    if (a[i] <= 0) phat <- phat + 1
      else phat <- phat
    i <- i + 1L
  }
  return(phat)
}

## basically what you had
## but instead of incrementally growing
## x, theta, and a (which requires R to make copies at each expansion
## instatiate them as n length vectors
buffon <- function(n) {
  x <- vector("numeric", n)
  theta <- vector("numeric", n)
  a <- vector("numeric", n)
  phat <- 0L
  i <- 1L
  while (i <= n) {
    x[i] <- runif(1, 0, 1/2)
    theta[i] <- runif(1, 0, pi/2)
    a[i] <- 1/2 * (sin(theta[i])) - x[i]
    if (a[i] <= 0) phat <- phat + 1
    i <- i + 1L
  }
  return(phat)
}

## take the trivial step of byte compiling the function
require(compiler)
cbuffon <- cmpfun(buffon)


## benchmark slow and uncompiled versus compiled
## on 500 iterations
set.seed(1)
system.time(for (i in 1:500) buffonslow(i))

set.seed(1)
system.time(for (i in 1:500) buffon(i))

set.seed(1)
system.time(for (i in 1:500) cbuffon(i))

## but maybe you will not be doing many iterations
## just one large number
set.seed(1)
system.time(buffonslow(50000))
set.seed(1)
system.time(buffon(50000))
set.seed(1)
system.time(cbuffon(50000))



On my system, here are the times.  Notice that just instatiating a
full sized vector yields a *huge* speedup when buffon is called for
large N.  Compiling helps further.
user  system elapsed
   9.82    0.00    9.96
user  system elapsed
   8.06    0.01    8.38
user  system elapsed
   3.35    0.00    3.44
user  system elapsed
  89.63    0.29   93.70
user  system elapsed
   3.21    0.00    3.31
user  system elapsed
   1.36    0.00    1.39

On Sat, Nov 19, 2011 at 9:08 PM, R. Michael Weylandt
<michael.weylandt at gmail.com> <michael.weylandt at gmail.com> wrote: