1

一些背景:我需要计算出从某个点到 3D 网格中每个单元格的距离,然后对该距离应用一个函数。我需要对多个点执行此操作,并在每个单元格中为所有点添加函数值。对于位于 (x,y,z) 的点,我可以使用以下代码执行此操作:

x <- c(1,2,3,4,5)
y <- x
z <- x
radius <- c(0.4,0.5,0.6,0.7,0.8)
numsphere <- length(x)
radius_buffer <- 0.2

xvox <- seq((min(x)-1),(max(x)+2),0.5)
yvox <- xvox
zvox <- xvox

probability_array <<- array(0,dim=c(length(xvox),length(yvox),length(zvox)))

for (j in 1:length(yvox)){ # for every y element
  for (i in 1:length(xvox)){ # for every x element
    for (k in length(zvox):1){ # for every z element
      for (n in 1:numsphere){ # for the total number of points
        dist_sd <- ((xvox[i]-x[n])^2+(yvox[j]-y[n])^2+(zvox[k]-z[n])^2)^0.5
        probability_array[i,j,k] <- probability_array[i,j,k] + 
                                    round(exp(-1*(dist_sd-radius[n])^2/(2*radius_buffer^2)),3)
          }
        }
      }
    }

输出是一个数组,绘制的结果如下所示:

probability_array <- probability_array/max(probability_array)
contour3d(probability_array,level=c(0.2,0.8,0.9),x=xvox,y=yvox,z=zvox,color = c("aquamarine","gold","darkorange"),alpha = c(0.1,0.2,0.5),add=T)

输出图形

我试图将其并行化,因为它似乎很适合它,但无法让它工作。我试过了:

cl<-makeCluster(detectCores(),type="SOCK")
registerDoSNOW(cl)

for (j in 1:length(yvox)){
  for (i in 1:length(xvox)){
    for(k in length(zvox):1){
      probability_array[i,j,k] <- foreach(n=1:numsphere, .combine='+') %dopar% {
        dist_sd <- ((xvox[i]-x[n])^2+(yvox[j]-y[n])^2+(zvox[k]-z[n])^2)^0.5
        round(exp(-1*(dist_sd-radius[n])^2/(2*radius_buffer^2)),3)
      }
    }
  }
}

和类似的东西:

r <- foreach(j=1:length(yvox)) %:% foreach(i=1:length(xvox)) %:% foreach(k=length(zvox):1) %:% foreach(n=1:numsphere, .combine='+') %do% {

        dist_sd <- ((xvox[i]-x[n])^2+(yvox[j]-y[n])^2+(zvox[k]-z[n])^2)^0.5
        probability_array[i,j,k] <- probability_array[i,j,k] + round(exp(-1*(dist_sd-radius[n])^2/(2*radius_buffer^2)),3)
        probability_array[i,j,k]

}

但我错过了一些重要的东西。任何帮助将不胜感激。干杯

4

1 回答 1

4

当并行化计算时,由于它引入的开销,最好并行运行大块计算,而不是小的 - 外部循环,而不是内部循环。

但是,在这种情况下,不需要并行化计算:您可以将它们向量化。

# 3-dimensional analogue of row() and col()
dim3 <- function( a, i ) { 
  stopifnot( length(dim(a)) == 3 )
  r <- a
  if( i == 1 ) { r[] <- rep(1:dim(a)[1], dim(a)[2] * dim(a)[3]) }
  if( i == 2 ) { r[] <- rep(1:dim(a)[2], each = dim(a)[1], times = dim(a)[3]) }
  if( i == 3 ) { r[] <- rep(1:dim(a)[3], each = dim(a)[1] * dim(a)[2]) }
  r
}

probability_array <- array(0,dim=c(length(xvox),length(yvox),length(zvox)))
i <- dim3(probability_array,1)
j <- dim3(probability_array,2)
k <- dim3(probability_array,3)
for (n in 1:numsphere){
  dist_sd <- sqrt(
    (xvox[i]-x[n])^2 + (yvox[j]-y[n])^2 + (zvox[k]-z[n])^2
  )
  probability_array <- probability_array + 
    # Rounding intermediate results looks suspicious
    round(exp(-1*(dist_sd-radius[n])^2/(2*radius_buffer^2)),3)
}
于 2013-08-21T08:22:24.010 回答