Реферат: Варианты алгоритма возведения в степень повышение точности и ускорение

Как, никто этого еще не придумал?

Не берусь судить. Вероятно, задача о том, как максимально быстро возвести действительное положительное число в произвольную действительную степень, решалась примерно столь же часто, как и вставала, - а вставала, полагаю, не раз. И все же не так давно я с ужасом обнаружил, что RTL из состава Borland Delphi последних версий (как Delphi 6, так и Delphi 7) подходит к решению не более профессионально, чем прилежный пятиклассник, впервые столкнувшийся с такой проблемой.

Взглянем на исходный код функции Power из модуля Math, любезно предоставленный Borland Software:

function Power(const Base, Exponent: Extended): Extended; begin if Exponent = 0.0 then Result := 1.0 { n**0 = 1 } else if (Base = 0.0) and (Exponent > 0.0) then Result := 0.0 { 0**n = 0, n > 0 } else if (Frac(Exponent) = 0.0) and (Abs(Exponent) <= MaxInt) then Result := IntPower(Base, Integer(Trunc(Exponent))) else Result := Exp(Exponent * Ln(Base)) end;

Примечательно, что в благих целях оптимизации процессор оставляют наедине с целой толпой ветвлений, приводящих его, в конце концов, в общем случае к пресловутому решению пятиклассника, а именно, к тривиальной формуле

(1) x**y = exp(ln(x**y)) = exp(y*ln(x)).

Здесь x**y означает возведение x в степень y, a exp(x) = e**x.

Что плохого в таком подходе к решению? Во-первых, в набор инструкций FPU не входит ни операция вычисления exp(x), ни взятия натурального логарифма ln(x). Соответственно, результат вычисляется в несколько этапов, в то время как можно пойти более прямым путем, как будет показано ниже. За счет этого падает скорость вычисления; кроме того, здесь действует интуитивное правило, которое грубо можно сформулировать так: чем больше операций выполняется над числом с плавающей запятой в регистрах сопроцессора, тем больше будет и суммарная погрешность результата.

ПРИМЕЧАНИЕ Позднейшая проверка показала, что как Visual C из Visual Studio .NET, так и C++ Builder 4.5 реализуют возведение в степень более качественно. Используемый в них вариант концептуально не отличается от того решения, которое я хочу предложить.

Есть предложение

Давайте проведем инвентаризацию. Какие инструкции сопроцессора связаны с возведением в степень или взятием логарифма? Приведу краткую выдержку из [1] и [2]:

F2XM1 – вычисляет 2**x-1, где -1<x<1.

FSCALE (масштабирование степенью двойки) - фактически считает 2**trunc(x), где trunc(x) означает округление к 0, т.е. положительные числа округляются в меньшую сторону, отрицательные – в большую.

FXTRACT – извлекает мантиссу и экспоненту действительного числа.

FYL2X – вычисляет y*log2(x), где log2(x) – логарифм x по основанию 2.

FYL2XP1 – вычисляет y*log2(x+1) для -(1-1/sqrt(2))<x<1-1/sqrt(2) c большей точностью, нежели FYL2X. Здесь sqrt(x) означает вычисление квадратного корня.

Вот, в общем-то, и все. Но уже на первый взгляд этого хватает, чтобы понять, что задача может быть решена более прямо, чем предлагает RTL Borland Delphi.

Действительно, почему не заменить показатель степени в (1) на 2? Ведь неперово число отнюдь не является родным для двоичной арифметики! Тогда получится

(2) x**y = 2**log2(x**y) = 2**(y*log2(x)).

--> ЧИТАТЬ ПОЛНОСТЬЮ <--