Я написал программу, которая работает с задачей 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
}
}
}
Спасибо!
Поскольку вам нужно перебирать значения x, вы не можете векторизовать это. В какой-то момент R должен работать с каждым значением x отдельно и по очереди. Возможно, вы сможете запустить вычисления на отдельных ядрах ЦП, чтобы ускорить процесс, возможно, используя foreach в одноименном пакете.
В противном случае (и это просто скрывает цикл от вас) оберните основную часть вашего цикла как функцию, например:
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.
Я вывернул это «наизнанку», создав вектор x, где i-й элемент является значением после каждой итерации алгоритма. Результат относительно понятен, т.к.
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
}
Это можно оптимизировать, чтобы (а) отслеживать только те значения, которые все еще находятся в игре (через idx) и (б) избегать ненужных операций, например, ifelse оценивает оба аргумента для всех значений x, x/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, а затем выполнить безоговорочное деление на два.
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
Кажется, что на этапе деления на два можно сделать более крупные шаги, например, 8 — это 2^3, поэтому мы можем сделать 3 шага (прибавить 3, чтобы считать) и закончить, 20 — это 2^2 * 5, поэтому мы можем сделать два шага и перейти к следующей итерации в 5. Реализации?
Другой подход признает, что часто приходится пересматривать низкие числа, так почему бы не запомнить их и сэкономить на пересчете?
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-кратным ускорением? Кроме того, рекурсия может стать слишком глубокой, но рекурсия может быть записана как цикл (что, вероятно, в любом случае быстрее).
