Итак, зададимся целью вычислить 1000000! (безо всяких приближений) за вменяемое время (порядка 1-2 минут). Предполагается, что работа с длинными числами встроена в язык, причём их умножение реализовано более-менее адекватно (ну хотя бы алгоритмом Карацубы). Далее в посте будет использоваться Ruby, в котором именно Карацуба и реализован.
Попытки посчитать "по определению" ни к чему толковому не приводят:
def factorial(n)
r = 1
for i in 1..n do
r *= i
end
return r
end
user system total real 25000! 3.320000 0.060000 3.380000 ( 3.393819) 50000! 14.210000 0.190000 14.400000 ( 14.458802) 100000! 58.970000 0.870000 59.840000 ( 60.115690)Как видим, время ~ n2, хотя число множителей ~ n. Невесело как-то. Эдак 1000000! нам часа два ждать придётся...
А дело в том, что почти все операции — это умножение длинного числа на короткое, в то время как алгоритмы перемножения длинных чисел ведут себя лучше всего при перемножении чисел равной длины. Доказывать это теоретически мне лень, проще это на практике показать.
Чтобы сделать длины перемножаемых чисел примерно равными, поступим следующим образом. Заметим, что
Здесь оба сомножителя, как нетрудно понять, примерно равной длины. С каждым из них можно проделать тот же трюк. Получаем следующую рекурсивную реализацию:
def g(n, step)
return n if n <= step
return g(n, step*2) * g(n-step, step*2)
end
def factorial(n)
g(n, 1)
end
Ничего сложного в этом нет, а время выполнения улучшилось на порядок: user system total real 25000! 0.330000 0.000000 0.330000 ( 0.365312) 50000! 1.030000 0.030000 1.060000 ( 1.098816) 100000! 3.560000 0.030000 3.590000 ( 3.666322)И в принципе, если лень заморачиваться, это вполне неплохой алгоритм. Тот самый 1000000! вычисляется за 4 с небольшим минуты.
Этот алгоритм можно ещё немножко ускорить, если заметить, что
class Factorial
@shift = 0
def g(n, step)
return n if n <= step
if n.even? and step == 2 then
@shift += n/2
return g(n/2, 1)
end
return g(n, step*2) * g(n-step, step*2)
end
def fac(n)
@shift = 0
v = g(n, 1)
return v << @shift
end
end
Однако можно ускорить вычисление факториала ещё раза эдак в два! Но для этого нужно применять куда более хитрые методы.(На этом вступление заканчивается ^____^)
Для начала вспомним, как возводить число в степень:
def power(a, p)
return 1 if p.zero?
t = power(a, p/2) ** 2
t *= a if p.odd?
return t
end
Применив совсем нехитрые соображения, мы уменьшили число умножений от O(n) до O(log n).На похожей идее основан самый быстрый на сегодняшний день алгоритм вычисления факториала под названием PrimeSwing. Однако никакого описания алгоритма мне нагуглить не удалось (потому-то я и пишу всю эту муть, хе-хе).
Так вот, ключевая идея: будем вычислять n! как произведение двух сомножителей.
Первый:
Второй — соответственно всё остальное, т.е.
(U/L — целое, так как число сочетаний из 2(n div 2) по (n div 2) целое.)
Так что сосредоточимся на вычислении второго сомножителя. Предположим, что у нас имеется таблица простых чисел, меньших n (до миллиона она вычисляется за доли секунды). Рассмотрим все простые числа
Случай 1:
Случай 2:
Случай 3:
Случай 4:
Вот теперь код по вышеприведённой ссылке должен быть более-менее понятен.
Алгоритм, написанный на Ruby, вычислил 1000000! за 84 секунды. (UPD от 30.05.2011: в 1.9.3-dev появился Toom-Cook, и теперь время составляет 47 секунд. Круть!) Для сравнения: GMP, одна из самых скоростных библиотек для работы с числами такой величины, вычисляет это число за 7 секунд. Однако там, в отличие от Ruby, используются куда более хитрые методы для умножения длинных чисел, например, алгоритм Шёнхаге-Штрассена. Так что результат весьма неплох.
UPD от 31.07.2011: в связи с интересом общественности к данному посту решил выложить и свой код на Ruby — ещё один язык исполнения лишним не будет =)
https://gist.github.com/1116938
No comments:
Post a Comment