АЛГОРИТМ 706

Интерполяция по Эйткену [Е1]

Если даны (n+l) абсцисс Хо, Xi,..., Хп и (п+1) значений функции f{xo), f{xi),..., }{Хп), то процедура aitken строит полином Лагранжа ff{x) степени я такой, что ff{Xi)-f{xi). Следовательно, для любого данного значения хх получается значение ff{xx). Для построения полинома ff{xx) применяется итеративная схема Эйткена. Поскольку массив f используется для хранения промежуточных значений, то его начальное значение затирается. Процедура aitken пригодна как для равных, так и для неравных интервалов [Хг,Хц.1].

procedure aitken(x,f,n,xx) result: (ff);

value xx,n; real xx,ff; integer n; array x,f; begin integer i,j;

for j:=0 step 1 until n-1 do for i:=]-}-l step 1 until n do

f[i]: = ( (XX-4J]) X f[il- (xx-x[i]) X f[j]) / (хт-хШ) \ ff: = f[n] end aitken;

Свидетельство к алгоритму 706

Алгоритм 706 является стереотипным переизданием алгоритма 70а, который успещно используется сейчас на многих мащинах, и в частности на БЭСМ-6.

При пользовании этим алгоритмом не следует задавать больших значений п, это может привести к переполнению или потере точности. Обычно достаточны значения п=2 или п-3. Поэтому обычно приходится дополнять обращения к такой процедуре операторами, обеспечивающими выбор из данной таблицы трех или четырех точек таких, чтобы точка XX оказалась в интервале [xo,Xm].

Свидетельство к алгоритму 70а

Алгоритм 70а получен в результате исправления, сокращения и ординарной переработки алгоритма 70 (Mifsud С. J. «САСМ», 1961, № II).

Алгоритм 70а проверен с помощью транслятора ТА-1 для функции е*. При хх=й.ЪЪ для равных интервалов х получено значение ff= 1.78605373 (по справочнику [1] f/= 1.7860). При xx=Q.S& для не-)авных интервалов х получено ff=2.36311017 (по справочнику [1] 7=2.3632).

Алгоритм 696 не пябликуется здесь, потому что соответствующий алгоритм 69 «Прослеживание цепочки» (Maych И. «САСМ», 1961, № 9) не был подтвержден расчетами ни в журнале «САСМ», ни авторами данного выпуска.

Подтверждение к алгоритму 70

А. Релф (Relph А. Р. «САСМ», 1962, № 7)

Алгоритм 70 был транслирован с использованием транслятора DEUCE ALGOL и дал удовлетворительные результаты после добавления ... *

Сзидетельстзо к алгоритму 716 [G6]

Алгоритм 716 не публикуется здесь, потому что соответствующий алгоритм 71 («САСМ», 1961, № И) содержал ошибки и был в дальнейшем заменен в журнале «САСМ» алгоритмами 86, 87, 102, 130, 202, дающими перестановки при -более экономном расходовании машинной памяти и машинного времени.

АЛГОРИТМ 726

Генератор композиций [А1]

Процедура сотрХ (сокращение от composition - композиция) преобразует данную композицию с из й частей целого положительного числа п ;в следующую по порядку композицию, если таковая существует. Если же процедура сотр\ выполняется с каждой новой -композицией после ее получения, то образуются все композиции при условии, что в качестве первоначальной была взята композиция 1, 1,.. .,1, я- Если с имеет форму я-й-Ы, 1, I,..., 1, то следующей композиции не существует и с заменяется вектором из к нулей [3, с 129].

procedure compl (с,к);

value к; integer к; integer array с; begin integer j;

j:=k; test: if j>l then begin ifc[j]=lthen

begin j: = j-1; go to test end; cO-l]: = c[j-l]+l; c[k]: = clj]-l; if jkthen c[j]: = l end else

for j: = l step 1 until к do c[j]: = 0 end comp 1;

Процедура comp2 образует и выдает на печать по каналу k всевозможные композиции из р частей целого положительного числа п. procedure comp2(n,p,k);

value n,p,k; integer n,p,k; begin integer j; integer array c[l:p];

for j: = l step 1 until p-1 do c[j]: = l;

о[р]:=п-Р-Ы; start:outarray(k,c);

j:=p;

* Указываются две поправки и возможность сокращения алгоритма 70, учтенные в алгоритме 70а. (Прим. ред.)

test: if j>l then begin if c]j] = l then

begin j:=j-1; go to test end; c[j-l]: = c[j-l]+l; c[p]: = c[i]-l; if jp then clj]:==l; go to start end end comp2;

Свидетельство к алгоритму 726

Алгоритм 726 отличается от алгоритма 72а только тем, что в процедуре сотр2 оператор с меткой «start:» был заменен на

start: outarray(k,c);

Алгоритм 726 был транслирован в системе ТА-IM для М-220 с примером, приведенным в «Свидетельстве к алгоритму 72а», и дал такие же результаты.

Свидетельство к алгоритму 72а

Алгоритм 72а написан в двух вариантах заново в более краткой, экономичной и наглядной форме, чем алгоритм 72 (Не Herman L., Ogden S. «САСМ», 1961, № 11). Перевод «Подтверждения к алгоритму 72» (Collison D. М. «САСМ», 1962, № 8) здесь не приводится, поскольку оно потеряло свое значение после написания алгоритма 72а.

Алгоритм 72а проверен для п=5 и k=3. Получено шесть следующих трехкомпонентных композиций числа 5:

1112 2 3 с 12] 12 3 12 1 с 13} 3 2 1 2 1 1

АЛГОРИТМ 736

Неполные эллиптические интегралы [S21]

Процедура ellint вычисляет значения неполных эллиптических интегралов первого и второго рода F{phi,k) и E{phi,k), где phi берется в радианах. Если ii>l или \phi\>-nl2, то Е и F приравниваются 10 иначе они будут содержать вычисленные интегралы. Алгоритм использует константы п/2= 1.57... и tg{l) = 1.557393 ...

Более полную информацию об этом алгоритме можно получить в работе Дидонато и Хершея [21].

procedure ellint (k,phi) result: (ee,ff);

value k,phi; real k,phi,ee,ff; begin real kp,sn,cs,hl,h2,al,a2,sl,s2,s3,s4,pl,p2,ml,m2,nl,n2,

tl,t2,dl,d2,d3,d4; integer n;

sl: = s2:=s3: = s4:=0; n:=0;

hl: = l; sn: = sin(phi); cs: = cos(phi);

if abs(kxsn) < 1.557393 then go to small;