s0 <- 1
m <- .07
r <- .03
s <- .2
dt <- .0001
I <- list()


me <- (r - .5 * s^2)*dt
sd <- s * sqrt(dt)


for(i in 1:1000) {
  prices <- list()
  prices[length(prices)+1] <- s0
  T <- dt
  while (T<=2) {
    if(length(prices) == 2/dt) break
    w <- rnorm(1,0,1)
    w2 <- me + w * sd
    prices[length(prices)+1] <- as.numeric(prices[length(prices)]) * exp(w2)
    T <- T + dt
  }
  I[length(I)+1] <- mean(as.numeric(prices))
  print(i)
}

print("Price Estimation:")
print(exp(-1*r*dt) * mean(as.numeric(I)))

print("Error estimation:")
print(sqrt(var(as.numeric(I))/1000))