Индийский алгоритм возведения в степень

"Индийский алгоритм" возведения в степень

Этот алгоритм вычисления натуральной n-й (n > 0) степени целого числа x выглядит совсем просто:

при n = 1 xⁿ = x;
при n > 1 xⁿ = x^{n mod 2}((x^{n
div 2})².

Основная цель этого алгоритма - сократить количество умножений при возведении в степень. Например, по этому алгоритму x⁵ = x·(x²)², т.е. достаточно трех умножений вместо четырех: x·x·x·x·x. Одно умножение экономится за счет того, что x² хранится как промежуточное значение и умножается само на себя. Точно так же x¹⁰ = 1·(x⁵)² = (x⁵)², что требует лишь четырех умножений (три из них для вычисления x⁵) вместо девяти "лобовых". Но здесь придется хранить сначала x², а потом x⁵.
Очевидно, что вычисление xⁿ сводится к вычислению x^{n div 2}, запоминанию результата, возведению в квадрат и умножению его на x при нечетном n. Итак, вычисление xⁿ описывается рекурсивной функцией:

function pow(x, n: integer) : integer;
var t: integer;
begin
	if odd(n) then t:= x else t:= 1;
	if n = 1 then pow:= x
	else pow:= t*sqr(pow(x, n div 2))
end;

    Как видим, промежуточные сомножители хранятся в локальной памяти процессов выполнения вызовов функции, а именно, в переменных, поставленных в соответствие ее имени.
    Теперь опишем зависимость глубины рекурсии вызовов функции от значения аргумента. В каждом следующем вложенном вызове значение аргумента n меньше предыдущего значения, по крайней мере, вдвое. Поскольку при n = 1 происходит возвращение из вызова, то таких уменьшений значения аргумента n не может быть больше чем log₂n. Следовательно, глубина рекурсии вызова с аргументом n не превышает log₂n.
    Такую глубину можно считать хорошим свойством алгоритма. При каждом выполнении вызова происходит не более одного деления, возведения в квадрат и умножения, поэтому общее количество арифметических операций не больше 3log₂n. При больших значениях n это существенно меньше "лобовых" n - 1 умножений. Например, при n = 1000 это примерно 30.
    Отметим, что при некоторых значениях n приведенный алгоритм не дает наименьшего количества умножений, необходимых для вычисления n-й степени. Например, при n = 15 по этому алгоритму необходимо выполнить 6 умножений, хотя можно с помощью 3-х умножений вычислить x⁵, после чего помножить его на себя дважды (всего 5 умножений). Однако написать алгоритм, который задает вычисление произвольной степени с минимальным количеством умножений, - не совсем простая задача.
    Построим нерекурсивный аналог приведенного алгоритма. Представим вычисление по рекурсивному алгоритму в таком виде:
x¹³ = (x⁶)² · x¹ = ((x³)² · x⁰)² · x¹ = (((x¹)² · x¹)² · x⁰)² · x¹
    Этому соответствует следующая обработка вычисляемых показателей степеней:
13 = 6 · 2 + 1 = (3 · 2 + 0) · 2 + 1 = ((1 · 2 + 1) · 2 + 0) · 2 + 1
    Как видим, вычислению степеней соответствует вычисление значения 13, представленного полиномом относительно 2. Коэффициентами его являются цифры двоичного разложения числа 13. Нетрудно убедиться, что вычислению степени с произвольным показателем n точно так же соответствует вычисление числа n, представленного двоичным разложением. Причем это разложение-полином записано по схеме Горнера. Раскроем в нем скобки:
1 · 2³ + 1 · 2² + 0 · 2¹ + 1 · 2⁰
    Коэффициенты при 2⁰, 2¹, 2² и так далее - это последовательные остатки от деления на 2 чисел:
n, n div 2, (n div 2) div 2 и так далее.
    Причем остатку 1 соответствует в рекурсивном алгоритме присваивание t:= x, а 0 - присваивание t:= 1. Таким образом, двоичное разложение, например, числа 13 по степеням двойки соответствует такому представлению x¹³:
x^2³ · x^2² · 1 · x^2⁰
    Итак, достаточно вычислять степени:
x^2⁰ = x, x^2¹ = x², x^2² = (x²)², x^2³ = (x^2²)² и так далее
и соответствующие им остатки от деления на 2 показателей:
n, n div 2, (n div 2) div 2, ((n div 2) div 2) div 2 и так далее
накапливая в произведении лишь те степени двойки, которые соответствуют остаткам 1. В следующем алгоритме произведение степеней накапливается в переменной t, а степени двойки - в переменной x:

function pow(x, n: integer): integer;
var t: integer;
	notfin: boolean;
begin
	t:= 1;
	notfin:= true;
	while notfin do
	begin
		if odd(n) then t:= t*x;
		n:= n div 2;
		if n > 0 then x:= x*x else notfin:= false
	end;
	pow:= t;
end;

Сайт управляется системой uCoz