О суммировании ряда Лагерра методом Эйлера-Кноппа в задаче обращения преобразования Лапласа Текст научной статьи по специальности «Математика»

УДК 518:517.432.1

Вестник СПбГУ. Сер. 1, 2008, вып. 4

М. М. Кабардов

О СУММИРОВАНИИ РЯДА ЛАГЕРРА МЕТОДОМ ЭЙЛЕРА—КНОППА В ЗАДАЧЕ ОБРАЩЕНИЯ ПРЕОБРАЗОВАНИЯ ЛАПЛАСА

1. Введение

Методам обращения преобразования Лапласа посвящена обширная литература. Одно семейство методов состоит в разложении оригинала в ряд по специальным функциям, в частности по многочленам Лагерра. Этот метод был разработан Пиконе [1] и Трикоми [2] и получил развитие у многих авторов (см. [3-8]). В настоящее время имеются различные варианты реализации этой схемы в зависимости от критериев выбора произвольных параметров, приемов ускорения сходимости и т. д. Мы рассматриваем метод обращения преобразования Лапласа с представлением оригинала в виде ряда по многочленам Лагерра, который затем суммируется методом Эйлера—Кноппа с целью ускорения сходимости. В работах [6-8] предложены способы выбора параметра суммирования методами Эйлера—Кноппа в нескольких частных случаях расположения особенностей изображения Лапласа. При этом выбор параметра р суммирования ограничен условием р < 1, р € И.

В настоящей работе изучается задача выбора оптимального значения параметра р в общем случае комплексного числа.

Итак, пусть дано изображение

/•то

Е (в) = в-8'/(г) ¿г, Ие в>А (А € И). Jo

Предположим, что вЕ(в) регулярна в бесконечности. Оригинал разложим в ряд по многочленам Лагерра:

то

/(г) = ^2 (ьг), (1)

к=0

где Ь — положительный параметр, Ь > 2А. Коэффициенты ак определяются из тейлоровского разложения функции

к=0

Преобразование Эйлера—Кноппа последнего ряда имеет вид

то то ь к /и\

<р(г) = = 5>(Р) , Ак(р) = ]Г (-р)^а, .

к=0 к=0 ^ Р ' ]=0

Множество допустимых значений параметра р , при которых область сходимости преобразованного ряда содержит круг сходимости исходного (такое свойство преобразования ряда будем называть регулярностью), было подробно рассмотрено в [9]. Множество соответствующих значений параметра р обозначим через М.

© М.М. Кабардов, 2008

2. Преобразование Эйлера—Кноппа

Теорема. Пусть функция-оригинал f (t) представима в виде (1), а параметр p удовлетворяет неравенству Rep < 1. Тогда

ff.s ( bpt\^ Ak(p) r f bt \ Доказательство. Известно [10] представление полиномов Лагерра

„t fOO

Lk{t) = yj e~vVkJo (2Vvt) d'q . (3)

Под знаком интеграла стоит целая функция, и интеграл от нее вдоль произвольного луча, исходящего из начала координат под углом таким, что | arg п| < п/2, сходится к Lfc(i). В самом деле, известно [11] асимптотическое разложение Jo(z):

Jo w . (А)1/2 Ls -1) f,-!,.^ - (, -1) E(-d^)

v ' \ s=0 s=0 J

при z ^ ж в секторе | arg z| < n — S (S > 0). О т с юда

JoW=o(|.r1/2|cos(z-^)|)=o(N-V2eN) .

Тогда подынтегральное выражение в (3) есть величина порядка

что, очевидно, влечет сходимость интеграла вдоль луча с углом | arg п| < п/2.

Из этих рассуждений и леммы Жордана следует, что при замене п = т/(1 — p) интеграл (3) не меняется при дополнительном требовании Reт/(1 — p) >0 (или, что то же самое, Rep < 1). Поэтому

гоо / \ / / J

-т fc____( pT \ Т (о I Tt

Аргумент функции .То выбирается так, что Р1е >0. Подставляя в последний ин-

теграл разложение

, pt \ ^ 1 ( Pt

и используя выражение (3), получаем

Подставив это в (1) и изменив порядок суммирования, получим формулу (2).

Замечание. В доказательстве используется ограничение Rep < 1. Оно выполнено при условиях, наложенных на изображение и параметр суммирования. Действительно, если sF(s) регулярна в бесконечности и суммирование регулярно, то |p | < r < 1 (см. [9]).

то

Следуя [6], назовем коэффициентом сходимости (КС) произвольного ряда ^ gn(z)

n=0

величину R = lim sup |gn(z)|1/n и рассмотрим вопрос выбора параметра p так, чтобы

n—>то

КС ряда (2) был минимальным.

3. Выбор параметра суммирования

Прежде всего заметим, что КС ряда (2) совпадает с КС ряда

^ Afc(p)

k=0 (1 - p)fc+1'

Это следует из асимптотической формулы для многочленов Лагерра [10] при k ^ ж:

Lk(z) = ^=e^(-kz)-^exp{2(-kz)^) + О (*r«/2)j ,

где z принадлежит комплексной плоскости с разрезом вдоль положительной полуоси, Cj (z) регулярны в этой же плоскости с разрезом. Отсюда lim sup |Lk(z)|1/k = 1 при

k—>то

любом значении z (при вещественных z > 0 см. [6]), что и делает очевидным наше утверждение.

Пусть t0 = argmax |tj|, tj = Sj/(sj — b), Sj —особенности sF(s). j

Утверждение 1. Суммирование Эйлера—Кноппа с параметром p регулярно тогда и только тогда, когда выполнены условия:

1) p = at0, 0 < а < 1,

2) |p | + A(p) = |t0|, A(p)= lim |Ak (p)|1/k = max |j — p |.

k—>то j

Необходимость этих условий установлена в [9], а достаточность очевидна. Легко убедиться, что из второго условия вытекает первое. Но приведенная формулировка удобнее для использования в дальнейшем.

Утверждение 2. Справедливо равенство

■ A(p) . „. .

arg mm--- = arg mm А (p .

& м \l-p\ м yi'

Из утверждения 1 имеем

А(Р) = (1-оот ..

\1-р\ |1-о*°| ' [ '

Так как А(р) = (1 — а)|г0| достигает минимума на М при максимальном а € [0,1], не выводящем р = аг0 из М, достаточно показать, что (5) убывает при 0 < а < 1.

Заметим, что в силу требований на изображение и параметр суммирования неравенство |г0| < 1 выполняется, тем самым знаменатель (5) на [0,1] в нуль не обращается. Найдем производную по а и приравняем ее нулю:

o -|l-at°|2-(l-a)(a|¿0|2-Ret0) _

11 |l-aí0|5/2

Отсюда находим, что единственная стационарная точка а0 = (1 — Re t0)/(Re t0 — |t012) не принадлежит отрезку [0,1]. Из того, что производная ((1 — а)/| 1 — at0|) непрерывна, не обращается в нуль на [0,1] и ((1 — а)/|1 — at0|)' = |t0|(Ret0 — 1) < 0, получаем

a=0

требуемое.

Было установлено [9], что p и A(p) суть соответственно центр и радиус замкнутого круга K(p, A(p)), который содержит все особенности tj, причем по крайней мере одна особенность лежит на границе этого круга.

Отсюда следует, что A(p)/|1 — p | = sin ß/2, где ß — угол, под которым круг K(p, A(p)) виден из точки t = 1. Поэтому геометрически задача arg min A(p)/|1 — p | состоит в нахождении центра круга K(p, A(p)), который виден под наименьшим углом из точки t = 1. Из утверждения 2 следует, что для этого (при требовании регулярности) достаточно найти круг наименьшего радиуса, который целиком лежит в круге K(0, |t01) и содержит все точки tj (см. рис. 1).

Рис. 1. Случай регулярного р.

В частном случае двух особых точек эта задача имеет решение

_ |£°|2- \Ь\2

где особенности обозначены так, что |£° | > 1.

4. Решение задачи в исходной плоскости

При отображении 5 = Ь£/(£ — 1) две прямые, проходящие через точку £ = 1 под углом 7 друг к другу, перейдут в прямые, проходящие через точку 5 = Ь под углом —7. Внутренность окружности дК(р, А(р)) перейдет во внутренность ее образа. Прямые, проходящие через £ = 0, перейдут в окружности, проходящие через точки 5 = 0, 5 = Ь. Окружность |£| = |£° | перейдет в окружность с диаметром [Ь|£°|/(|£°| — 1), Ь|£°|/(|£°| + 1)] .

С учетом этих замечаний задача нахождения оптимального параметра р, обеспечивающего регулярное суммирование Эйлера—Кноппа, может быть решена в плоскости (5). Для этого достаточно выполнить следующие действия.

1. Найти одну особенность s функции sF(s), лежащую на границе круга B, обладающего следующими свойствами:

— B содержит все особенности sF(s) и имеет наименьший возможный диаметр;

— диаметром служит отрезок вида [bc/(c — 1), bc/(c + 1)], c G (0,1).

2. Найти круг Bi (so , r) наименьшего радиуса, обладающий следующими свойствами:

— круг B1 (so , r) касается окружности dB в точке S;

— B1 (so ,r) содержит все особенности sF(s).

3. По радиусу r и центру s0 = xo + iy0 найти s1 и s2:

b — xo — ixo , .

51,2 = So ±ra, а = —-:-. (6)

|b — xo — ixo |

4. Вычислить ропт = (t(s1) + t(s2))/2, t(s) = s/(s — b).

Эти действия отображены на рис. 1.

Если мы откажемся от регулярности и поставим задачу безусловной минимизации min A(p)/|1 — р то круг K(p, A(p)) не будет в общем случае принадлежать K(0, |to|). Но решение рОпт удовлетворяет требованию Re рОпт < 1 теоремы. В плоскости (s) эта задача решается следующим образом.

1. Найти круг, который содержит все особенности sF(s) и виден под наименьшим углом из точки s = b (см. рис.2).

Рис. 2. Случай нерегулярного p.

2. По радиусу r и центру s0 найти s1 и s2 по формуле (6).

3. Вычислить р'опт = (t(si) + t(s2))/2, t(s) = s/(s — b).

4. Численные эксперименты

Расчеты показали, что выбор параметра p по описанной схеме (в комплексной плоскости) дает существенный выигрыш в скорости сходимости по сравнению с выбором на вещественной оси. Ниже приведены результаты сравнения скорости сходимости и точности при различных значениях параметра суммирования (см. табл. 1).

Далее p' обозначает оптимальный параметр без условия регулярности, p'' — оптимальный параметр, обеспечивающий регулярность суммирования, p''' — оптимальный параметр на вещественной оси, p = 0 соответствует обычному суммированию. В полях таблицы стоят величины max |f (xj) — Sn(xj)|/|f (xj)|, xj = jT/100, j = 0,1,..., 100

j

(Sn —отрезки ряда (2)).

Коэффициенты Ак (р) вычисляются с использованием быстрого преобразования Фурье (БПФ) по следующей схеме. Строится интерполяционный многочлен для функции Ф(т) = 1/(1 + рад) (т/(1 + рт)) = ^Акр)тк. В качестве узлов интерполяции {wj} берутся узлы Вандермонда Wj = ехр (2пг/К^ — 1/2)), ] = 1, 2,..., N (они обеспечивают равномерную сходимость интерполяционного процесса в единичном круге). Для оптимизации числа операций в алгоритме БПФ выбрано N = п . Затем, в качестве Ак(р) берутся первые п коэффициентов БПФ.

Вычисления проведены в среде MATLAB с относительной погрешностью округлений £ « 10-16. Тестовый пример имеет следующие исходные данные: /(¿) = ви + ^(в) = 1/(в — г) + 1/в2, Ь =1, Т = 6.

Таблица 1.

р' р" р'" р = 0

р 1.46е—001—г3.53е—001 2.50е—001 - г2.50е—001 1.55е—015 0

А(р)/\1-р\ 4.1421е—001 4.4721е—001 7.0711е—001 7.0711е—001

п= 15 6.2269е—004 7.3171е—004 8.7007е—003 8.7007е—003

га=20 1.7364е—005 3.0281е—005 9.6129е—004 9.6129е—004

п=25 4.3314е—007 9.9969е—007 1.8374е—004 1.8374е—004

Литература

1. Picone M. Sulla transformazione di Laplace // Rend. Atti. Accad. Naz. Lincei. 1935. Vol. 21. P. 306-313.

2. Tricomi F. Transformazione di Laplace e polinomi di Laguerre // Rend. Atti. Accad. Naz. Lincei. 1935. Vol. 13. P. 232-239.

3. Shohat J. Laguerre polynomials and the Laplace transform // Duke Math. J. 1940. Vol. 6. P. 615-626.

4. Piessens R., Branders M. A. Numerical inversion of the Laplace transform using generalized Laguerre polynomials // Proc. IEE. 1971. Vol. 118. P. 1517-1522.

5. Davis B., Martin B. Numerical inversion of the Laplace transform: a survey and comparison of methods // J. Comput. Phys. 1979. Vol.33. P. 1-32.

6. Лебедева А. В., Рябов В. М. Об обращении преобразования Лапласа с помощью рядов Лагерра и квадратурных формул // Методы вычислений. Вып. 19. СПб., 2001, C. 123-139.

7. Gabutti B., Lepora P. The numerical performance of Tricomi's formula for inverting the Laplace transform // Numer. Math. 1987. Vol. 51. P. 369-380.

8. Gabutti B., Lyness J. N. Some generalizations of the Euler—Knopp transformation // Numer. Math. 1986. Vol. 48. P. 199-220.

9. Кабардов М. М. О применении метода суммирования Эйлера—Кноппа к ряду Лагерра // Методы вычислений. Вып. 22. СПб., 2008. C. 77-81.

10. Сегё Г. Ортогональные многочлены. М., 1962. 500 с.

11. Олвер Ф. Асимптотика и специальные функции. М., 1990. 528 с.

Статья поступила в редакцию 18 мая 2008 г.