2

私は 3n + 1 問題 (別名「不思議な数」など) で動作するプログラムを作成しました。しかし、それには二重のループがあります。どうすればベクトル化できますか?

コードは

count <- vector("numeric", 100000)
L <- length(count)

for (i in 1:L)
{
x <- i
   while (x > 1)
   {
   if (round(x/2) == x/2) 
     {
     x <- x/2
     count[i] <- count[i] + 1 
     } else
     {
     x <- 3*x + 1
     count[i] <- count[i] + 1
     }
   }
}

ありがとう！

4

3 に答える 3

9

i 番目の要素がアルゴリズムの各反復後の値であるベクトル x を作成することで、この「裏返し」を変えました。結果は比較的わかりやすい

f1 <- function(L) {
    x <- seq_len(L)
    count <- integer(L)
    while (any(i <- x > 1)) {
        count[i] <- count[i] + 1L
        x <- ifelse(round(x/2) == x/2, x / 2, 3 * x + 1) * i
    }
    count
}

これは、(a) まだ有効な値のみを (idx 経由で) 追跡し、(b) 不必要な操作を避けるように最適化できます。たとえば、ifelse は x のすべての値に対して両方の引数を評価し、x/2 は 2 回評価されます。

f2 <- function(L) {
    idx <- x <- seq_len(L)
    count <- integer(L)
    while (length(x)) {
        ix <- x > 1
        x <- x[ix]
        idx <- idx[ix]
        count[idx] <- count[idx] + 1L
        i <- as.logical(x %% 2)
        x[i] <- 3 * x[i] + 1
        i <- !i
        x[i] <- x[i] / 2
    }
    count
}

元の関数 f0 を使用すると、次のようになります。

> L <- 10000
> system.time(ans0 <- f0(L))
   user  system elapsed 
  7.785   0.000   7.812 
> system.time(ans1 <- f1(L))
   user  system elapsed 
  1.738   0.000   1.741 
> identical(ans0, ans1)
[1] TRUE
> system.time(ans2 <- f2(L))
   user  system elapsed 
  0.301   0.000   0.301 
> identical(ans1, ans2)
[1] TRUE

微調整は、奇数値を 3 * x[i] + 1 に更新してから、無条件に 2 で割ることです。

x[i] <- 3 * x[i] + 1
count[idx[i]] <- count[idx[i]] + 1L
x <- x / 2
count[idx] <- count[idx] + 1

これを f3 とすると (なぜ今朝 f2 が遅いのかわかりません!)

> system.time(ans2 <- f2(L))
   user  system elapsed 
   0.36    0.00    0.36 
> system.time(ans3 <- f3(L))
   user  system elapsed 
  0.201   0.003   0.206 
> identical(ans2, ans3)
[1] TRUE

2 で割る段階で、より大きなステップを実行できるようです。たとえば、8 は 2^3 であるため、3 ステップ (カウントに 3 を追加) を実行して終了できます。20 は 2^2 * 5 なので、次のことができます。 2 つのステップを実行し、5. 実装で次の反復に入ります。

于 2010-12-18T23:24:39.400 に答える
4

の値を反復処理する必要があるため、xこれを実際にベクトル化することはできません。ある時点で、R は x の各値を個別に順番に処理する必要があります。foreachおそらく同じ名前のパッケージを使用して、別の CPU コアで計算を実行して速度を上げることができる場合があります。

それ以外の場合 (これはループを隠しているだけです)、ループの本体を関数としてラップします。

wonderous <- function(n) {
    count <- 0
    while(n > 1) {
        if(isTRUE(all.equal(n %% 2, 0))) {
            n <- n / 2
        } else {
            n <- (3*n) + 1
        }
        count <- count + 1
    }
    return(count)
}

sapply()次に、一連の数値に対して関数を実行するために使用できます。

> sapply(1:50, wonderous)
 [1]   0   1   7   2   5   8  16   3  19   6  14   9   9  17  17
[16]   4  12  20  20   7   7  15  15  10  23  10 111  18  18  18
[31] 106   5  26  13  13  21  21  21  34   8 109   8  29  16  16
[46]  16 104  11  24  24

または、それ自体が関数であるVectorizeベクトル化されたバージョンを返すために使用することもできます。wonderous

> wonderousV <- Vectorize(wonderous)
> wonderousV(1:50)
 [1]   0   1   7   2   5   8  16   3  19   6  14   9   9  17  17
[16]   4  12  20  20   7   7  15  15  10  23  10 111  18  18  18
[31] 106   5  26  13  13  21  21  21  34   8 109   8  29  16  16
[46]  16 104  11  24  24

現時点で、標準の R ツールで得られるのはこれくらいだと思います。@Martin Morgan は、 R のベクトル化された機能を使用した問題を独創的に解決することで、これよりもはるかに優れたことができることを示しています。

于 2010-12-18T22:54:01.893 に答える
2

別のアプローチでは、低い数値を頻繁に再訪することを認識しているので、それらを覚えて再計算のコストを節約してみませんか？

memo_f <- function() {
    e <- new.env(parent=emptyenv())
    e[["1"]] <- 0L
    f <- function(x) {
        k <- as.character(x)
        if (!exists(k, envir=e))
            e[[k]] <- 1L + if (x %% 2) f(3L * x + 1L) else f(x / 2L)
        e[[k]]
    }
    f
}

これは

> L <- 100
> vals <- seq_len(L)
> system.time({ f <- memo_f(); memo1 <- sapply(vals, f) })
   user  system elapsed 
  0.018   0.000   0.019 
> system.time(won <- sapply(vals, wonderous))
   user  system elapsed 
  0.921   0.005   0.930 
> all.equal(memo1, won) ## integer vs. numeric
[1] TRUE

これはうまく並列化できないかもしれませんが、50倍のスピードアップでは必要ないのではないでしょうか。また、再帰が深くなりすぎる可能性がありますが、再帰はループとして記述される可能性があります（とにかく、おそらくより高速です）。

于 2010-12-20T03:02:11.090 に答える