泊松分布函数无限循环

Question

我有以下代码，试图模拟金属原子沉积到冷基板上;但是，它在一个无限循环中运行。

有人能看出我在哪里出错吗?

l <- 20
n <- 2000
e <- 1000
lsize <- matrix(0,l,l)
deposits <- rep(0,n)
avg.deposits <- rep(0,n)
prob <- rep(0,n)
n.deposits <- rep(0,n)
for(m in 1:e){
  for(j in 1:l){
    for(k in 1:l){
      lsize[j,k] <- 0
    }
  }
  for(i in 1:n){
    ra <- runif(1)
    x <- floor(1+l*ra)
    ra <- runif(1)
    y <- floor(1+l*ra)
    lsize[x,y] <- lsize[x,y]+1
    s <- 0
    for(j in 1:l){
      for(k in 1:l){
        if(lsize[j,k] <- 1){
          s <- s+1
        }
      }
    }
    n.deposits[i] <- n.deposits[i]+s
  }
}
for(i in 1:n){
  avg.deposits[i] <- n.deposits[i]/e
  prob[i] <- avg.deposits[i]/(l*l)
  deposits[i] <- i
}
plot(deposits, prob)

Answer 1

不存在无限循环问题。

如果你继续用较小的l, n, e参数运行你的代码，这很容易检查。当增加上述任何参数时，您的代码的扩展都不是最优的(在本例中是超线性的)。

明显的点:

1. Preallocate矩阵。每个循环不要重复分配lsize。
限制函数调用;到runif()。您不必多次调用同一个函数。在循环外调用它一次，生成你想要的随机数，然后在循环内访问该行的下一个元素。
2. 使用print和cat语句打印您使用的循环计数器。尝试小的值，以确保他们的程序做你想要的，然后将计数器设置为数千。
3. 尽可能对代码进行矢量化。如。如果'a = runif(100)'并且您希望设置a < 0.5等于4的所有实例，则没有理由依次遍历a的所有元素。a[ a < 0.5] = 4就够了。