LI. Инициализация

L4. Вывод Wn

Рис. 17. Алгоритм L обращения степенного ряда. Если алгоритм L применить к примеру z - t-t, можно вычислить следующее.

1 2 3 4 5

-1 О

Uo Ui

2 3 4

6 10 15

20 35

1 1 2

В упр. 8 показано, что небольшая модификация алгоритма L позволит решать значительно более общие задачи, прилагая только чуть больше усилий.

Рассмотрим сейчас решение уравнения

Uiz + U2Z + Usz + ... = t + V2t + V3t + (14)

относительно t и получим коэффициенты степенного ряда

t = Wiz + W2Z+ W3Z+ W4Z*+ . (15)

Уравнение (10)-это частный случай (14) при Ui = 1, U2 = U3 = = О- Если Ui ф О, можно предположить, что U\ - 1, заменив z на ((/iz), но мы рассмотрим общее уравнение (14), поскольку U\ может равняться нулю.

Алгоритм Т {Обобщенное, обращение степенных рядов). В этом интерактивном алгоритме вводятся значения Un и F„ из (14) и выводится значение Wn из (15) для п = 1, 2, 3, ..., N. При вычислениях используется вспомогательная матрица Ттп, l<m<n<N.

Т1. [Инициализация.] Присвоить п <- 1. Занесем два первых введенных значения (а именно - (7i и Fi) в Гц и 14 соответственно. (Должно выполняться равенство Vl = 1.)

Т2. [Вывод Wn.] Вывести значение Tin, которое равно И„.

ТЗ. [Ввод Un, Vn.] Увеличить п на 1. Если п > N, алгоритм заканчивает работу, иначе - запоминает два следующих введенных значения (а именно - Un и Vn) в Тщ и Vn.

Т4. [Умножение.] Присвоить

Ттп TiiTm-l.n-l +Ti2Tm-l,n-2 Н-----Ь Tin-m+lTm-l,m-l

и Tin Tin - VmT,nn ДЛЯ 2 < m <п. (После этого шага для 1 < m < п получим

t" = Ттшг + Г„,ш+1г"+1 + • • • + TmnZ + 0(z"+i). (16)

Легко проверить (16) индукцией по m > 2 и, когда т - 1, получить Un = Tin + У2Т2П Н-----Ь VnTnn согласно (14) и (16).) Возвратиться к шагу Т2.

Соотношение (16) объясняет механизм этого алгоритма, предложенного Генри К. Тэчером (мл.) (Henry С. Thacher, Jr.) [САСМ 9 (1966), 10-11]. Время счета алгоритма, по существу, такое же, как и у алгоритма L, но требуется значительно больший объем памяти для хранения данных. Пример работы алгоритма приведен в упр. 9.

Другой подход к обращению степенных рядов предложен в работе R. Р. Brent and Н. Т. Kung, JACM 25 (1978), 581-595. Он основан на том факте, что стандартные итерационные процедуры, которые используются для нахождения корней уравнений с действительными числами, можно также применять к уравнениям для степенных рядов. В частности, можно рассмотреть метод Ньютона для приближенного вычисления действительного числа t, такого, что f{t) = О, а заданная функция / хорошо ведет себя в окрестности t: если х является хорошим приближением к t, то ф{х) = X - f{x)/f{x) будет даже лучшим приближением. Записав х = t + е, получим f{x) = fit) + ef{t) + 0{е), f{x) = f{t) + 0(e); следовательно, ф{х) = t + е - [О + ef{t) + 0{e))/{f{t) + 0(e)) = t + 0{e). Применим эту идею к степенным рядам. Пусть f{x) = V{x) - U{z), где U и V - степенные ряды из (14). Найдем степенной ряд t от z, такой, что f{t) = 0.

Пусть X = Wiz Н-----h Wn-iz"~ =t + 0(г") - "приближение" к t порядка п, тогда

ф{х) = X - f{x)/f{x) будет приближением порядка 2п, поскольку для этих / и t выполняются предположения метода Ньютона.

Другими словами, можно воспользоваться следующей процедурой.

Алгоритм N {Обобщенное обращение степенного ряда методом Ньютона). Данный "полиинтерактивный" алгоритм вводит значения Un а Vn согласно (14) для 2* < п < 2*+ и выводит значения Wn согласно (15) для 2* < п < 2*+, получая ответы группами по 2* значений одновременно для к = О, 1, 2, ..., К.

N1. [Инициализация.] Присвоить iV 1. (Получим N = 2.) Ввести первые коэффициенты Ui и Vi {Vi = 1) и присвоить Wi <- Ui.

N2. [Вывод.] Вывести Wn для N <п < 2N.

N3. [Ввод.] Присвоить N <- 2N. Если N > 2, алгоритм закончил работу, иначе - ввести значения Un и У„ для N <п < 2N.

N4. [Шаг Ньютона.] Воспользуемся алгоритмом для композиции степенных рядов (см. упр. 11), чтобы вычислить коэффициенты Qj и Rj {О < j < N) степенного ряда

UiZ + --- + U2N-lZ - V{Wiz + + Wn iz-)

= Rz"" + Riz + + Rm-iz"- + 0(2), V{Wiz + • • • + Wn-iz"-) = Qo + Qiz + --- + Qn-iz- + 0{z),

где V{x) = i + Vba; Н----и V{x) = I + 2V2XA----. Затем возьмем Wj,W2N-1

в качестве коэффициентов степенного ряда

и возвратимся к шагу N2.

Время работы данного алгоритма при получении N = 2 коэффициентов равно T{N), где

T{2N) = T{N) + (время на выполнение шага N4) + 0(N). (17)

Прямые алгоритмы для композиции и деления на шаге N4 имеют порядка шагов, значит, алгоритм N работает медленнее алгоритма Т. Тем не менее Брент (Brent) и Кунг (Kung) нашли метод, которым требуемая композиция степенных рядов выполняется с помощью 0(iV log iV)/ арифметических операций (в упр. 6 приведен алгоритм для деления, работающий еще быстрее). Таким образом, в (17) показано, что обращение степенных рядов можно выполнить только с помощью 0{N log iV)/ операций, когда iV -> оо. (С другой стороны, константа пропорциональности такова, что N должно быть действительно большим, прежде чем алгоритмы L и Т перестанут относиться к "быстродействующим" методам.)

Историческая справка. Ж. Н. Брамхел (J. N. Bramhall) и М. А. Чаппл (М. А. Chappie) впервые опубликовали метод обращения степенных рядов, требующий 0{N) операций, в САСМ 4 (1961), 317-318,503. Это, по существу, автономный алгоритм, эквивалентный методу, который приведен в упр. 16, с таким же временем счета, как у алгоритмов L и Т.

Итерация рядов. Чтобы изучить поведение итеративного процесса i„ f{xn-i), следует изучить п-кратную композицию данной функции / саму с собой, т. е. Хп = /(/(... f{xo) ...))• Определим /М(а;) = х п /W(i) = /(/t"-4(a;)) так, что

/t"+"l(a;) = /Н(/И(х)) (18)

для всех целых т, п > 0. Во многих случаях обозначение /"(х) имеет смысл и для отрицательных целых п. Если /["1 и - взаимно обратные функции, а именно - X = /["](/t~"l(2;)), и если обратные функции определены единственным образом, то (18) справедливо для всех целых тпп. Обращение рядов - это, по существу, операция нахождения обратного степенного ряда /[""(i). Например, соотношения (10) и (И) устанавливают, что z - V[W{z)) и t = W{y{t)); таким образом, W =

Предположим, что заданы два степенных ряда, V{z) = z + Уг + • • и W{z) -

z + Wz Л----, таких, что W = Пусть и - любая не равная нулю постоянная.

Рассмотрим функцию

V{z) = W{uV{z)). (19)

Легко видеть, что (/([/(г)) = И(п2у(г)) и вообще

C/I"](z) = И(цЧ/(г)) (20)

для всех целых п. Следовательно, имеем простое выражение для п-й итерации (/"5, которую можно вычислить приблизительно с одинаковыми затратами труда