CC BY
26
УДК 519.65 Вестник СПбГУ. Сер. 10. 2016. Вып. 4
А. Ю. Утешев, И. И. Боровой
РЕШЕНИЕ ЗАДАЧИ РАЦИОНАЛЬНОЙ ИНТЕРПОЛЯЦИИ С ИСПОЛЬЗОВАНИЕМ ГАНКЕЛЕВЫХ ПОЛИНОМОВ
Санкт-Петербургский государственный университет, Российская Федерация, 199034, Санкт-Петербург, Университетская наб., 7—9
Работа посвящена задаче построения рационального интерполянта
r(x) = p(x)/q(x), {r(xj ) = yj , {xj ,yj }N=1 CC, {p(x), q(x)}C С [x] .
В развитие результата К. Якоби интерполянт представляется в виде отношения ганкеле-вых полиномов, т. е. полиномов вида Нк(x) = det[ci+j—1 — Ci+j—2x]Kj=1. Порождающая последовательность {ck}kg^ выбирается в виде {£]N=1 xkyj/W'(xj)}fcgN для полинома
q(x) и N=1 xk/(yj W'(xj))}кеы для полинома p(x); здесь W(x) = П N=1(x — xj). Приводятся условия разрешимости задачи и несократимости получаемой дроби. В дополнение к формальному построению решения в детерминантной форме в настоящей статье предложена процедура эффективного вычисления соответствующих ганкелевых полиномов. Она основана на тождестве Якоби—Йоахимшталя, связывающем ганкелевы полиномы трех последовательных порядков линейным соотношением вида
аНк(x) — (x + в)Нк—1 (x) + 1/аНк—2 (x) = 0
при некоторых константах {а, в} C С. Доказательство этого соотношения также приводится в статье вместе с дополнительным обсуждением вырожденного случая а = 0. На основании изложенных результатов может быть развернута процедура вычисления ганкелевых полиномов, рекурсивная по их порядку. Такая возможность позволяет получить не только интерполянт с фиксированными степенями полиномов p(x) и q(x), но и все семейство интерполянтов при различных комбинациях степеней: deg p + degq ^ N — 1. Библиогр. 12 назв.
Ключевые слова: рациональная интерполяция, ганкелевы матрицы и полиномы, алгоритм Берлекампа—Месси.
A. Yu. Uteshev, I.I. Baravy
SOLUTION OF THE RATIONAL INTERPOLATION PROBLEM VIA THE HANKEL POLYNOMIAL CONSTRUCTION
St. Petersburg State University, 7—9, Universitetskaya nab., St. Petersburg, 199034, Russian Federation
The problem of rational interpolant construction is treated as
r(x) = p(x)/q(x), {r(x j ) = yj }N=1, {xj ,yj CC, {p(x), q(x)}C С [x] .
At the basis of one result by C. Jacobi, the interpolant is represented as a ratio of two Hankel polynomials, i. e. polynomials of the form Нк(x) = det[ci+j—1 — Ci+j—2x]Kj=1. The generating sequence for these polynomials is selected as {£]N=1 xkyj/W'(xj)}fcgN for q(x) and as {2N=1 xj/(yjW'(xj))}fcgN for polynomial p(x); here W(x) = J^N=1 (x — xj). The conditions
Утешев Алексей Юрьевич — доктор физико-математических наук, профессор; alexeiuteshev@ gmail.com
Боровой Иван Иванович — аспирант; [email protected]
Uteshev Alexei Yur'evich — doctor of physical and mathematical sciences, professor; alexeiuteshev@ gmail.com
Baravy Ivan Ivanovich — postgraduate student; [email protected] © Санкт-Петербургский государственный университет, 2016
for the solubility of the problem and irreducibility of the obtained fraction are also presented. In addition to formal representation of the solution in determinantal form, the present paper is focused also at the effective computational algorithm for the Hankel polynomials. It is based on a little known identity by Jacobi and Joachimsthal connecting a triple of the Hankel polynomials of successive orders:
аНк(x) — (x + @)Hk-i(x) + 1/аНк-2 (x) = 0 ,
here {а, в} С С are some constants. The proof of this relation is also contained in the paper along with discussion of a degenerate case а = 0. With these results, a procedure for the Hankel polynomial computation can be developed which is recursive in its order. This gives an opportunity not only to compute a single interpolant with specialized degrees for p(x) and q(x) but also to compose the whole set of interpolants for an arbitrary combination for the degrees: deg p + deg q ^ N — 1. The results of the paper can be applied for problems of Approximation Theory, Control Theory (transfer function reconstruction from frequency responses) and for error-correcting coding (Berlekamp—Welch algorithm). Although the presented results are formulated for the case of infinite fields, they are applicable for finite fields as well. Refs 12.
Keywords: rational interpolation, Hankel matrices and polynomials, Berlekamp—Massey algorithm.
1. Введение. Рассмотрим задачу интерполяции для рациональных функций. Задача. Для таблицы значений переменных x и y
I N
X Х\ Х2 xN
У У1 У2 Ум
{хз }]=1 все различны, (1)
построить ее рациональный интерполянт, т. е. рациональную функцию такую, что
{Фз ) = у0 }?=1 . (2)
Здесь г(х) = р(х)/д(х) при р(х) = рохп + р1хп-1 + ... + рп, д(х) = дохт + д1хт-1 + ... + Чт, Ро = 0, до = 0, и N = п + т +1.
Замечание 1.1. В дальнейшем не будем делать различия между интерпо-лянтами, числитель и знаменатель которых домножены на общий числовой множитель.
Исторически первое решение задачи было предложено Коши в 1821 г. — в развитие метода Лагранжа построения полиномиального интерполянта. Интерес к задаче возродился во второй половине ХХ в. и был связан прежде всего с приложениями в теории аппроксимации, теории управления (восстановление передаточной функции по частотным характеристикам) и помехоустойчивом кодировании: в последнем случае задача ставится в конечных полях [1]. В работах [2-4] обсуждаются как теоретические аспекты рациональной интерполяции, так и приложения к задачам вычислительной математики.
Следует упомянуть о некоторых особенностях задачи рациональной интерполяции, отличающей ее от интерполяции полиномиальной. В то время как последняя всегда имеет решение, задача интерполяции рациональной не всегда разрешима в указанной постановке. Для иллюстрации этого преобразуем (2) в систему уравнений
{Р(хз )= Уз Ч(хз )}^=1, (3)
или, в развернутом виде,
{pn +-----+ pixnn-1 + pox] = gmyj + qm-ixj yj +-----+ qixji-1yj + gox™ yj}N=1 , (4)
которая линейна относительно N +1 коэффициентов р(х) и д(х). Конструктивная разрешимость этой системы может быть установлена средствами линейной алгебры.
Пример 1.1. При degр(х) = 1, deg д(х) = 3 найти рациональный интерполянт для таблицы
X -1 0 1 2 3
У 1 1 1/3 3 1/13
Решение. Решая систему (4), получаем р(х) = х — 2, д(х) = х3 — х2 — х — 2. Однако р(2) = 0 и д(2) = 0, а потому условие г(2) = 3 не выполнено. Оно не будет выполнено, даже если сократить р(х) и д(х) на их общий линейный делитель. □ Природа этого феномена кроется в неэквивалентности перехода от (2) к (3), поскольку для некоторого узла хз может существовать решение линейной системы (4), удовлетворяющее обоим условиям: р(хз) =0 и д(хз) = 0. Вместе с тем в случае существования решение задачи может оказаться не единственным. Пример 1.2. Для таблицы
X -1 0 1 2 3
У 1 1 1/3 1/7 1/13
значений функции 1/(х2 + х + 1) существует бесконечно .много интерполянтов при degр(х) = 1, deg д(х) =3 в форме (х — Х)/((х — Х)(х2 + х + 1)) при X ^ { — 1, 0,1, 2, 3}. □
Настоящая статья посвящена подходу к решению задачи, основанному на идее, предложенной в 1846 г. Карлом Якоби [6]. Такой подход заключается в представлении полиномов числителя и знаменателя в виде определителей специального вида — так называемых ганкелевых полиномов. Элементами этих определителей, помимо мономов по х, являются некоторые рациональные симметрические функции элементов таблицы (1). Работа Якоби может считаться практически забытой в XX в.; ссылки на нее немногочисленны (и не всегда корректны). Возможно, данное обстоятельство связано с сомнениями в практической значимости предложенного подхода. Действительно, определитель большого порядка, зависящий от параметра, крайне неудобен для практических расчетов*. Однако же удалось обнаружить процедуру, позволяющую обойти упомянутую проблему: для ганкелевых полиномов существует рекурсивная по порядку определителей процедура их вычисления, сводящая расчет такого полинома к-го порядка к получению двух полиномов меньших порядков. Удивительно то, что автором идеи такой процедуры оказался тот же Якоби [6], а ее исчерпывающее обоснование было осуществлено его учеником Фердинандом Йоахимшталем [7]. Приведем этот результат в п. 2. В п. 3 излагается основной теоретический результат по разрешимости задачи рациональной интерполяции в терминах ганкелевых полиномов и на примере иллюстрируется применение результата Якоби—Йоахимшталя для получения всего семейства решений рассматриваемой задачи — при всех возможных комбинациях степеней числителя и знаменателя.
Замечание 1.2. В этой статье приводятся доказательства только основных результатов в сокращенном виде. Полные их версии изложены в [8].
2. Ганкелевы определители и полиномы. Для числовой последовательности (конечной или бесконечной)
{сз 3=о = {С0,С1,...} (5)
* Достаточно вспомнить проблему вычисления характеристического полинома матрицы.
матрица вида
1,3=1
со
С1
С1
С2
С2 С3
Ск-2 с— Ск Ск-1 Ск Ск+1
Ск-1 Ск
С2к-3 С2к-2
(6)
кк
называется ганкелевой матрицей порядка к, порожденной последовательностью (5). Ее определитель будем обозначать Нк({с}), или просто Нк, если это не вызывает путаницы.
Если заменим последнюю строку ганкелевой матрицы порядка к +1 строкой, составленной из степеней х, то соответствующий определитель
Нк(х; {с}) =
Со С1
С1 С2 С2 Сз
= ("1)А
Ск—1 Ск Ск+1
1 х 2 х2
С1 - Сох С2 - С1 х
С2 - С1Х сз - С2 х
Ск
Ск+1
С2к-1 к
х
(7)
(к+1)х(к+1)
Ск - Ск-1Х Ск+1 - СкХ
Ск - Ск-1Х Ск+1 - СкХ
С2к-1 - С2к-2Х
кк
или просто нк (х), называется к-м ганкелевым полиномом [9], порожденным последовательностью (5). Разложение определителя (7) по его последней строке дает
нк(х) = нкохк + нк1хк 1 + нк2хк 2 + ... при к
ко
Нк .
Таким образом, degНк(х) = к тогда и только тогда, когда Нк = 0.
Теорема 2.1 (Якоби, Йоахимшталь). Любые три последовательных ганкелевых полинома Нк—2(х), Нк—1 (х), Нк(х) связаны соотношением
Н2кНк-2 (х) + (Нккк—1,1 - Нк-1кк1 - НкНк-1х) Нк-1 (х) + Н—Нк(х) = 0. (8)
к
Доказательство*. Рассмотрим сначала случай, когда порождающая последовательность (5) задается в виде
т } 2к—1
Сз-^ 3
(9)
1=1 ) 3=0
при т > к и произвольных различных А1,..., Ат. Докажем справедливость следующих соотношений:
Т.А1Нк(Ае) = \ Н 6СДи • €{0,...,к - 1}, (Ю)
1 I Нк+1, если • = к. '
* Приводимое доказательство представляет собой переписанное в современных обозначениях оригинальное доказательство из [7].
Действительно,
J24 Hk А ) = £ Aj
со ci
ci C2
Ck-1 Ck
1
Aj
Ck
Ck+1 C2k-i
Ak
со ci
ci c2
ck i ck
cj
cj+1
ck
ck+1
c2k-1 cj+k
При з < к последний определитель обращается в нуль, поскольку содержит две одинаковые строки. При з = к он совпадает с Нк+1.
Пусть Нк-1 = 0 (т.е. degНк-1(х) = к — 1). Разделим Нк(х) на Нк-1 (х):
Нк(х) = Q(x)Hк-l(x) + Е(х). (11)
Здесь коэффициенты частного выражаются через коэффициенты Нк(х) и Нк-1(х):
гм \ /-> I Нк „ Нк-1Ък1 — Нккк-1,1 , ,
4{х) = Чо + при <31 = 77-, <?о =-779--• (12)
H
k-1
H 2
nk-1
Для нахождения коэффициентов остатка R(x) = Ro + R1X + • • • + Rk-2Xk 2 подставим x = A1,...,x = Am в (11):
{Hk (А,) = (Qo + QA) Hk-1(Ae) + ( Ro + RA + ••• + Rk-2Ahl-2) Ц
1
(13)
Суммирование этих равенств дает
т т т
^Нк(Хе) = QоJ2Нк-1(Хе) + ХеНк-1(Хе) + (едйо + + ••• + с-й-). 1
j=1
Согласно (10), получим
0 = coRo + c1R1 +-----+ ck-2Rk-2
Далее умножим каждое равенство (13) на соответствующие Х^ и просуммируем полученные равенства по I. При 3 € {1,...,к — 3} приходим к равенствам
0 = С3йо + Сз+1 Й1 + ••• + Сз+к-2Йк-2 .
При з = к — 2 имеем равенство слегка иного вида:
0 = НкQl + Ск-2Йо + Ск-1Й1 + ••• + С2к-4Йк-2 .
Объединяя полученные равенства в систему, рассматриваем ее как линейную относительно йо,..., йк-2 и разрешаем по формулам Крамера. Результат может быть записан как тождество
й(х)Нк-1 + HкQlHк-2 (х) = 0 ,
которое вместе с выражением (12) для Ql подтверждает истинность (8) для частного случая порождающей последовательности, заданной (9).
Рассмотрим теперь случай произвольной порождающей последовательности (5). Для любой заданной последовательности комплексных чисел С1,...,С2к-1 возможно
отыскать комплексные числа Ах ,...,А^ такие, что при I > 2к — 1 уравнения (9) будут совместными. Эти числа могут быть найдены как корни полинома степени £, первые 2 к — 1 сумм Ньютона [10] которого совпадают с {с^ }2=-1.
Для завершения доказательства следует заполнить пробел в аргументации предыдущего абзаца. В то время как числа сх,..., с2к—1 могут быть выбраны произвольными, число со оказывается целым положительным, а именно равным I. Таким образом, истинность (8) доказана лишь для со € N. Однако доказываемое равенство является алгебраическим относительно {с^ }2=о1 и, будучи выполненным для бесконечного набора целых чисел, должно выполняться для любого со € С. □ Тождество (8) позволяет организовать рекурсивную процедуру вычисления ган-келевых полиномов. В самом деле, предположим, что канонические представления полиномов Нк—2 (х) и Нк—1 (х) уже найдены:
Нк-1(х) = Нк—1,охк—1 + Нк—1,1хк—2 +-----+ Нк—1,к—1 при Нк—1,0 = Нк—1 .
В таком случае в (8) известны значения почти всех констант, за исключением Нк и Нк 1. Для последних имеются лишь детерминантные представления
Нк
к
со с1
с1 с2
ск—2 ск—1
ск—1 ск
ск—2 ск—1 ск 1 ск
с2к—4 с2к—3 с2к-3 с2к-2
Нк1 = —
со с1
с1 с2
ск — 2 ск—1 ск 1 ск
ск—2 ск
ск — 1 ск+1
с2к—4 с2к — 2
с2к-3 с2к-1
Эти определители отличаются от транспонированного детерминантного представления Нк—1(х) только последними столбцами. Разложения по элементам последних столбцов имеют одинаковые значения для соответствующих алгебраических дополнений, и потому формулы
Нко = Нк = ск—1Нк—1,к—1 + ск Нк—1,к—2 + ••• + с2к—2Нк—1,о, Нк1 = —(скНк — 1,к — 1 + ск+1Нк — 1,к — 2 + ••• + с2к—1 Нк — 1,о)
(14)
позволяют выразить Нко и Нк1 посредством уже известных коэффициентов полинома Нк—1(х).
Однако предложенный алгоритм рекурсивного вычисления Нк(х) не работает в случае, когда Нк—1 = 0. Модификация процедуры может быть осуществлена с помощью следующего результата.
Теорема 2.2. Пусть Нк—2 = 0, Нк—1 = 0. Если Нк—1,1 = 0, то
Н
Нк- 1(ж)=0 и Нк(х) = '"к2 Нк-2(х) ■
Нк 2
В противном случае
Нк-\{х) = ^ 1'1Нк-2(х)
Нк
к—2
(15)
Нк (х) =
Нк—2 0 0 Нк
НкНк—2Нк—1,1Нк—з(х) — Нк—2,1 Нк—2,2 2 х2 Нк—2 Нк—2,1 х 0 Нк—2 1 Нк1 Нк2 0 Нк — 2(х)
Нк-2
(16)
и
Замечание 2.1. Формулы теоремы 2.2 — (15) и (16) — позволяют производить рекурсивное вычисление Нк (х), когда полиномы Нк-2 (х) и Нк-з (х) уже посчитаны. Участвующие в этих формулах константы, такие как Нк-2,1, Ьк-2,2, Ьк-1,1 и Нк,1, либо считаются известными как коэффициенты ганкелевых полиномов, либо могут быть вычислены при помощи формул (14). Единственным исключением является Нк2. Для вычисления этой величины предлагается использовать формулу
Ьк2 —
(с2к-2 Ьк-2,0 + С2к-зЬк-2,1 + • • • + Си Нк-2,к-2У
и,
к-2
которая справедлива при Ик-1 = 0, Ик-2 = 0.
Замечание 2.2. Фактически, формулу (8) следует воспринимать в качестве первоисточника алгоритма, известного в настоящее время как алгоритм Берлекампа—Месси [11]; он был предложен для декодирования кодов Боуза— Чоудхури—Хоквингема и Рида—Соломона, а также для нахождения минимального полинома линейной рекуррентной последовательности.
3. Рациональная интерполяция. Основной теоретический результат статьи предварим следующей леммой [10]:
Лемма 3.1. Обозначим Ш(х) = (х — хз)• Справедливы равенства Эйлера— Лагранжа
N
Е
3=1
х)
Ш '(хз)
0, если к е{1,...,М — 2},
1, если к = N — 1.
(17)
Теорема 3.1. Пусть {уз = 0}^=г Вычислим последовательности
N
Тк
3=1
Ш '(хз)
}3=1
2т
к=0
Тк
К 1 хк
Е1 хз
3 = 1
Уз Ш'(хз)
2п-2
к=0
(18)
и построим соответствующие им ганкелевы полиномы Нт(х; {т}) и Нп(х; {т})^ Ес-
Ип({Т}) = 0
) N
(19)
(20)
{Нт(хз; {Т}) = 0}\=1
то существует единственное решение задачи рациональной интерполяции при degр(х) = п, deg д(х) ^ т = N — п — 1 Это решение можно представить в виде
р(х) = Ит+1({т})Нп(х; {т}) =
то
Т1
Тт — 1
Т1 Т2
Тт+1
Тт Тт+1
Т2т-1 Т2т
Т0 Т1
Тп-1
Т1 Т2
Тп
Тп+1 Т2п-1
(21)
и
ли
и
т
п
т
д(х) = Ип({Т})Нт(х; {т}) =
то
Т1
Т1 Т2
Тп— 1 Тп
Тп-1
Т2п-2
то
Т1
Т1 Т2
Тт-1 Тт
1х
Тт Тт+1
Т2т-1
хт
Доказательство. Если решение задачи существует, то равенства (4) справедливы. Домножим 3-е равенство на х,/Ш'(хз) при к € {0,. ..,т — 1} и просуммируем полученные равенства по 3. На основании равенств (17) приходим к системе уравнений
{Чт Тк + Чт-1Тк+1 + • • • + Ч1Тк+т-1 + Ч0Тк+т = 0}т=0 . Таким образом, знаменатель дроби должен удовлетворять соотношению
Ад(х) =
то
Т1
Тт — 1
Т1 Т2
Тт Тт+1
Т2т-1
хт
т
1
х
при некоторой константе А.
Подобным образом, умножая равенства (4) на х*/(узШ'(хз)) при I € {0,...,п — 1} и суммируя по 3, получим представление числителя в виде
Вр(х) =
то
Т1
1
Т1 Т2
Тп— 1 Тп
х
Тп
Тп+1 Т2п-1
= Ип({Т})хп +
при некоторой константе В. Согласно (19), имеем В = 0 и degр(х) = п.
Для нахождения множителей А и В подставим полученные выражения в (3):
{АНп(хз; {Т}) = ВузНт(хз; {т})}
N
з=1
(23)
Согласно (20), А = 0 и {Нп(хз; {Т}) = 0}^=1. Домножим каждое из равенств (23) на хт/Ш'(хз) и просуммируем получившиеся результаты. В силу свойства линейности определителя и с использованием (17) имеем
N £
Нп(
хз;{т})х
Ш '(хз
то
Т1
Тп-1 N
Т1 Т2
N хт+1 \ "
Ш'(хз) Ш'(хз)
Тп-1 Тп
Т2п-2
N хт+п-1
Х_А_
£
Ш '(хз)
Тп
Тп+1
Т2п-1 N хт+п
£
Ш '(хз)
п
т
х
з
Тп-1 Тп 0 0
Аналогичным образом получаем
то Т1
Т1 Т2
Тп-1
Тп
Тп
Тп+1
Т2п-2 Т2п-1 0 1
= Нп({Т}).
Е
3=1
Ж '(хз)
= Нт+1({г}).
Следовательно,
АНп({Т}) = БНт+1({г}).
Поскольку А = 0 и Нп({т}) = 0, то Нт+1({г}) = 0, и последнее равенство завершает доказательство представимости рационального интерполянта по формулам (21) и (22).
Доказательство того факта, что эти полиномы действительно обеспечивают выполнение равенств (3), более громоздко и основано на равенствах
N
Нт+1({Т }) = (-1Г V-1)/2Нп({Т})\{ Уз,
3 = 1
N
N
Нт(хк; {т}) = (-1^^-1)/2Нп(хк; {Т}) Ц
У3
) к=1
Замечание 3.1. Формулировка теоремы 3.1 принадлежит авторам статьи. Якоби в [5] предложил разыскивать знаменатель рациональной интерполянты (потенциального кандидата) в форме ганкелева полинома Нт(х; {т}) и после его вычисления свести задачу к случаю интерполяции полиномиальной. Он не рассматривал вопросы существования и единственности решения, а также не интересовался вычислительными аспектами практической реализации предложенного им подхода, включая процедуру, используемую при решении следующего примера и основанную на его же собственном результате, изложенном в п. 2.
Пример 3.1. Найти все рациональные интерполянты г(х) = р(х)/д(х), degр(х)+ deg д(х) ^ 6 для таблицы
X -2 -1 0 1 2 3 4
У 26/51 2 -1/2 1/6 -4/7 16/31 7/36
Решение. Поскольку максимально возможная степень числителей и знаменателей искомых интерполянтов равна 6, вычислим последовательности (18) для значений индексов {0,1,..., 12}:
то = -
897683 19123776
, ...,Т12
5257205447 2390472 '
то = -
2973 11648'
,Т12
1294306589 11648
Определим ганкелевы полиномы первого и второго порядков:
897683 119579 208609 2 321193 7649
Агцж; {г}) — ————— х+———, /г2(ж; {г}) — ———— х -
19123776 4780944 ^^^^ 50996736 50996736 1416576
"V" V V
/¿2,0 ^2,1 ^2,2
Расчет Нз(х; {т}) произведем с применением тождества (8):
Из(х; {г}) = - (V2 П1(х; {г}) + (х - + Н2(ж; М),
\П2,0 ] "2,0 \ "2,0 "з,0 /
где все константы, кроме "з,о = Из({т}) и Нз,1, уже известны. Для нахождения последних воспользуемся равенствами (14)
Нз,о = Из({т}) = таН2}о + Тз"2,1 + Т2Н22 = —4037/16998912 ,
Нз,1 = —(т5"2,о + Т4Н21 + тз"2,2) = 36263/50996736 .
Поэтому
4037 з 36263 2 767 41
7тч(ж; |тИ =--х Ч--х--х--.
^ '1 16998912 50996736 12749184 75888
Продолжая рекурсивное использование тождества (8), получим
п, . , 1915 4 1915 з 9575 1915
Тгл(х',\т\) =--х -\--х -\--х-\--,
^ 1 50996736 25498368 152990208 38247552'
, г 1Ч 1915 5 36385 4 6229 з 40711 2 359 3037 Н5(х]{т}) =--ж Н--ж Н--хл--х1--х+-
21855744 76495104 8999424 10927872 6374592 1195236
, г 1Ч 991 6 8887 5 3475 4 51575 з 8450 2 4892 416
Н6(х: |г|) =-хь--ж Н--ж Н--хл--х1--х-\--
^ 1 796824 1195236 2390472 1195236 298809 99603 42687
/<6,0 /<6,1 /<6,6
Необходимо выполнить еще одну итерацию, а именно вычислить
И7({Т }) = Т12Н60 + Т11 Н6,1 + •••+ Т6Н6,6 = —208/42687.
Таким образом, все возможные знаменатели интерполяционной дроби получены. Аналогичная рекурсивная процедура может быть организована для определения числителей:
2973 3037 1915 2 21065 1915
Нл (х: 1т|) =--ж -I--, Ло(х: 1т|) = - ж--ж +
Н1 { } — Н2 { } —
11648 11648' ' " ¿06496 745472, 372736'
Ь2,0=Н2({т}) /12,1 /¡<2,2
Нз,,0 = т1н2,О + ТзН21 + Т2Н22 = 5745/745472 , Нз,1 = —(Т5Н2,0 + Т4Н2,1 + ТзН2}2) = —28725/745472 ,
{г}) . - (Ь>V Н1(х; {г}) + Ы (х - Ы + Ы П2(ж; {7}) = \Н2,0) "2,0 \ "20 "з,0)
5745 3 28725 2 33771 369
_ Qß __I _ ,
745472 745472 372736 6656
, r~iN 36333 4 72771 3 11139 2 1005843 206523
ггл(х',\т\) =-x--x--x -\--x--,
v 745472 372736 28672 745472 372736'
, 625827 5 1708605 4 362367 3 2007361 2 3068941 119579 Hdx; ir}) =--xb H--x4--xö--xz -\--x -\--,
1 ,г " 745472 372736 745472 93184 186368 46592 '
Нб(х; {f}) =
897683 6 5805465 5 373613 4 24907053 3 9008491 2 392865 42687 — __ __^ _|___ _^ _ __
93184 93184 7168 93184 23296 23296 416
Наконец, составим рациональные интерполянты комбинациями найденных числителей и знаменателей:
, ^ Н7({т}) _
го, б(х —
Нб(х; {т}) 11648
2973 х6 - 17774 х5 + 3475 х4 + 103150 х3 - 67600 х2 - 117408 х + 23296 '
h6,oHi (x; {f}) _
ri,5 (х) =
hi,oH (х; {т}) 64(2973х - 3037)
13405 х5 - 72770 х4 - 105893 х3 + 569954 х2 + 8616 х - 388736 h5,oH2 (х; {f}) 7 х2 - 11 х + 2
Г24 I = ^- = -;--- ; • • • ;
,К> h2fiHA(x-{T}) Зж4—6ж3—5ж — 4 h!,oH6(х; {f}) _
r6,o (х
he,(
897683 6 1935155 5 4856969 4 8302351 3 9008491 2 130955 1 —___|_______|__^ _|____
19123776 6374592 19123776 6374592 4780944 1593648 2 '
□
В заключение коснемся вопроса, связанного с единственностью решения интерполяционной задачи. Какое из условий теоремы 3.1 блокирует ситуацию из примера 1.2, т. е. возможность существования общего нетривиального делителя числителя и знаменателя?
Следствие. Имеет место следующее равенство:
Нп({Т})Нт ({т }) = Нп+1({Т})Нт+1({т}).
Оказывается, что все четыре определителя из этого следствия напрямую связаны с объектом, известным как результант полиномов, т. е. с такой функцией коэффициентов этих полиномов, обращение которой в нуль гарантирует существование общего их корня [8].
4. Заключение. В настоящей статье развит основанный на идеях К. Якоби подход к решению задачи рациональной интерполяции, заключающийся в представлении решения посредством подходящих ганкелевых полиномов. Предложена процедура эффективного вычисления этих полиномов, позволяющая построить не только одиночную интерполянту, но и целое семейство интерполянт при всех возможных комбинациях степеней числителя и знаменателя. Этот аспект обеспечивает возможность выбора решения, не только формально удовлетворяющего интерполяционной
таблице, но и обладающего дополнительными свойствами, существенными для задач аппроксимации (как, к примеру, ограниченность на определенном интервале вещественной оси).
Следует отметить, что появление матриц ганкелевой структуры в решениях задач аппроксимации не должно считаться абсолютно неожиданным: достаточно вспомнить вид системы нормальных уравнений при построении полиномиальной аппроксимации интерполяционной таблицы по методу наименьших квадратов [10]. В книгах [9, 12] можно найти и другие задачи, в которых такие матрицы возникают естественным образом, — например, задачу интерполяции таблицы (1) комбинацией экспонент:
f (*о=£ m=ix.
Развитие предложенного подхода видится в направлении интерполяции многомерной. Кроме того, предполагается произвести сравнение вычислительной эффективности алгоритма с альтернативными подходами, в частности с представлением рационального интерполянта в барицентрической форме [3].
Литература
1. Berlekamp E., Welch L. Error correction of algebraic block codes. Patent USA 4633470. 1986.
2. Beckermann B., Labahn G. Fraction-free computation of matrix rational interpolants and matrix GCD's // SIAM J. Matrix anal. appl. 2000. Vol. 22, N 1. P. 114-144.
3. Berrut J.-P., Baltensperger R., Mittelmann H. D. Recent developments in barycentric rational interpolation // Intern. Series of numerical mathematics. Basel: Birkhauser, 2005. Vol. 151. P. 27-51.
4. von zur Gathen J., Gerhard J. Modern computer algebra. 2nd ed. Cambridge: Cambridge University Press, 2003. 785 p.
5. Jacobi C. G. J. Uber die Darstellung einer Reihe gegebner Werthe durch eine gebrochne rationale Function // J. reine angew. Math. 1846. Vol. 30. P. 127-156.
6. Jacobi C. G. J. De eliminatione variabilis e duabus aequationibus algebraicis //J. reine angew. Math. 1836. Vol. 15. P. 101-124.
7. Joachimsthal F. Bemerkungen aber den Sturm'schen Satz //J. reine angew. Math. 1854. Vol. 48. P. 386-416.
8. Uteshev A. Yu., Baravy I. Solution of interpolation problems via the Hankel polynomial construction. arXiv: cs.SC/1603.08752. 2016. 56 p.
9. Henrici P. Applied and computational complex analysis. Vol. 1. Power series, integration, conformal mapping, location of zeros. New York; London: Wiley and Sons, Inc., 1974. 700 p.
10. Фаддеев Д. К. Лекции по алгебре. М.: Наука, 1984. 416 с.
11. Блейхут Р. E. Быстрые алгоритмы цифровой обработки сигналов. / пер. с англ. И. И. Груш-ко. М.: Мир, 1989. 448 с. (Blahut R. E. Fast algorithms for signal processing.)
12. Иохвидов И. С. Ганкелевы и теплицевы матрицы и формы. М.: Наука, 1974. 263 с. (Iohvidov I. S. Hankel and Toeplitz matrices and forms: Algebraic theory.)
Для цитирования: Утешев А. Ю., Боровой И. И. Решение задачи рациональной интерполяции с использованием ганкелевых полиномов // Вестник Санкт-Петербургского университета. Сер. 10. Прикладная математика. Информатика. Процессы управления. 2016. Вып. 4. С. 31-43. DOI: 10.21638/11701/ spbu10.2016.403
References
1. Berlekamp E., Welch L. Error correction of algebraic block codes. Patent US 4633470, 1986.
2. Beckermann B., Labahn G. Fraction-free computation of matrix rational interpolants and matrix GCD's. SIAM J. Matrix Anal. Appl. 2000, vol. 22, no. 1, pp. 114-144.
3. Berrut J.-P., Baltensperger R., Mittelmann H. D. Recent developments in barycentric rational interpolation. Intern. Series of numerical mathematics. Basel, Birkhaauser Press, 2005, vol. 151, pp. 27-51.
4. von zur Gathen J., Gerhard J. Modern Computer Algebra. 2nd ed. Cambridge, Cambridge University Press, 2013, 785 p.
5. Jacobi C. G. J. Uber die Darstellung einer Reihe gegebner Werthe durch eine gebrochne rationale Function [On the representation of the set of given values via the rational function]. J. reine angew. Math. [Journal of Pure and Applied Mathematics], 1846, vol. 30, pp. 127-156.
6. Jacobi C. G. J. De eliminatione variabilis e duabus aequationibus algebraicis [On the elimination of
variable from two algebraic equations]. J. reine angew. Math. [Journal of Pure and Applied Mathematics], 1836, vol. 15, pp. 101-124.
7. Joachimsthal F. Bemerkungen iiber den Sturm'schen Satz [A note of Sturm's theorem]. J. reine angew. Math. [Journal of Pure and Applied Mathematics], 1854, vol. 48, pp. 386-416.
8. Uteshev A. Yu., Baravy I. Solution of interpolation problems via the Hankel polynomial construction. arXiv: cs.SC/1603.08752, 2016, 56 p.
9. Henrici P. Applied and computational complex analysis. Vol. 1. Power series, integration, conformal mapping, location of zeros. New York, London, Wiley and Sons Inc. Publ., 1974, 700 p.
10. Faddeev D. K. Lekcii po algebre [Lectures in Algebra]. Moscow, Nauka Publ., 1984, 416 p. (In Russian).
11. Blahut R. E. Fast algorithms for digital signal processing. New York, Cambridge University Press, 1985, 485 p. (Russ. ed.: Blahut R. Bystrye algoritmy cifrovoj obrabotki signalov. Moscow, Mir Publ., 1989, 448 p.)
12. Iohvidov I. S. Hankel and Toeplitz matrices and forms: algebraic theory. Boston, Birkhauser Press, 1972, 260 p. (Russ. ed.: Iohvidov I. S. Gankelevy i teplicevy matricy i formy. Moscow, Nauka Publ., 1974, 263 p.)
For citation: Uteshev A.Yu., Baravy I.I. Solution of rational interpolation problem via the Hankel polynomial construction. Vestnik of Saint Petersburg University. Series 10. Applied mathematics. Computer science. Control processes, 2016, issue 4, pp. 31-43.DOI: 10.21638/11701/ spbu10.2016.403
Статья рекомендована к печати проф. Л. А. Петросяном. Статья поступила в редакцию 30 июня 2016 г. Статья принята к печати 29 сентября 2016 г.
