Об обращении преобразования Лапласа некоторых специальных функций Текст научной статьи по специальности «Математика»

ОБ ОБРАЩЕНИИ ПРЕОБРАЗОВАНИЯ ЛАПЛАСА НЕКОТОРЫХ СПЕЦИАЛЬНЫХ ФУНКЦИЙ*

H. И. Порошина1, В. М. Рябов2

I. С.-Петербургский государственный университет, аспирантка, [email protected]

2. С.-Петербургский государственный университет, д-р физ.-мат. наук, профессор, [email protected]

1. Введение. Теория линейной наследственности (вязкоупругости) находит приложения в различных разделах математики, физики и в механике деформируемого твердого тела. В настоящее время вязкоупругие материалы применяются в различных отраслях техники. Отметим, что процесс создания элементов конструкций из них неотделим от процесса создания материала, из которого они изготавливаются. Поэтому, прежде всего, в силу дороговизны натурного эксперимента, важную роль играют математические модели названных процессов. В задачах линейной вязкоупругости [1], описывающих напряженное состояние на основе определяющего соотношения Больцмана—Вольтерра (пространственные координаты ниже для простоты опущены)

e(t) = ^ (a(t) +ßj0 K(t-T)a(T)dT^ , (1)

сравнительно просто находятся изображения по Лапласу решений (деформаций e(t) или напряжений <r(t)) с учетом соотношения (1). Сами решения медленно изменяются на конечном по t отрезке времени и допускают хорошие приближения вида ts-1Q(ta) при 0 < a < 1, s > 0, где Q(t) —некоторый многочлен. Такая ситуация характерна для длительных процессов деформирования [2]. Изображения таких функций равны p-sR(p-a), где R(x) —многочлен той же степени, что и Q(x).

Важнейшей задачей становится выбор подходящего ядра, описывающего поведение материала, и определение его параметров. Ядро K может иметь интегрируемую особенность в точке t = 0 [1]. Чаще всего в качестве такового берут дробно-экспоненциальную функцию Ю. Н. Работнова [1]

к=0

(ßt1+a)k Г((1 + ск)(1 + к))

3a(ß,<)=î“Ef7ïï^7^, -К а О. (2)

Способ определения параметров дробно-экспоненциальной функции по измеренной функции ползучести описан в работе [2].

Интеграл от этого ядра по полуоси £ ^ 0 должен быть конечным, для чего необходимо в < 0. Не умаляя общности, далее считаем в = -1, и пусть символ Эа(£) означает

Эа(-М).

В наследственной механике твердого тела наряду с функцией (2) широко используется и интеграл от нее с переменным верхним пределом. Для облегчения использования

* Работа выполнена при финансовой поддержке РФФИ (грант №08-01-00285). © Н. И. Порошина, В. М. Рябов, 2009

этих величин составлены таблицы функций [1, 3]

Г1(а,х) = £ аЭа(х), Г2(а,х) = £ а 1 / Эа(т) ¿т, х = 1а+1.

■! 0

Однако при решении конкретных задач необходимо вводить в память вычислительной машины части этих таблиц, соответствующие найденным параметрам Эа-функций, которые к тому же заранее неизвестны и определяются в процессе решения задачи (и в итоге таковых в таблице может не оказаться). При изменении параметров приходится эту работу проделывать заново, что неудобно и сопряжено с внесением ошибок.

Изображения по Лапласу функции Эа(£) и интеграла от нее ^ Эа(т) ¿т равны, соответственно,

1 1

------, —------а = 1 + а.

ра + 1 р(ра + 1)

По этим изображениям определяются образы напряжений а(Ь) и деформаций е(Ь) из уравнения (1).

Естественным обобщением дробно-экспоненциальных функций и ядер Гаврильяка— Негами [2] являются ядра, изображения которых имеют вид

МаЛ7 (А,р)= С[(ра+1 + в)7 + А]-1. (3)

В частных случаях отсюда получаем целую экспоненту (а = А = 0, 7 = 1), ядро Абеля (в = А = 0, 7 =1), ядро Работнова (А = 0, 7 = 1), ядро Гаврильяка—Негами (А = 0).

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

Поскольку составление таблиц таких ядер, аналогичных таблицам Эа-функций, нецелесообразно ввиду большого количества входных параметров, практическое использование этих ядер для аппроксимации экспериментальных данных и при решении краевых задач вязкоупругости возможно лишь при разработке достаточно точных и эффективных методов вычисления этих функций и определяемых ими искомых решений.

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

В данной работе рассматриваются методы обращения преобразования Лапласа в предположении, что заданное изображение Г (р) искомой функции-оригинала фактически зависит от ра, где а — произвольное положительное число из интервала (0,1) (в частности, изображение может иметь имеет вид (3)). В случае а =1 получаются известные методы [4], в противном случае получаем новые формулы, обладающие большей точностью по сравнению с известными для определенного класса изображений и имеющие большое прикладное значение.

2. Деформация контура интегрирования. Задача обращения интегрального преобразования Лапласа состоит в нахождении решения уравнения

Г(р) = / е рх/(х) ¿х,

0

в котором F(p) —известное изображение, f (x) —искомый оригинал. Для простоты будем считать, что функция F(p) регулярна в полуплоскости Re (p) > 0, чего всегда можно добиться домножением оригинала на соответствующую экспоненту. Как правило, точное обращение осуществить не удается, и потому возникает необходимость разработки и применения приближенных методов. Наиболее полно возможные подходы к задаче обращения и их реализация описаны в книге [4]. В первую очередь, следует назвать построение квадратурных формул для приближенного вычисления интеграла Римана—Меллина

1 rc+itt

№ = L-\F)(t) = — eptF(p) dp, с > 0, (4)

J c-itt

задающего обращение преобразования Лапласа.

Напомним, что интеграл (4) понимается в смысле главного значения, он не зависит от с ив случае разрыва оригинала в точке t мы получаем полусумму предельных значений оригинала слева и справа от точки t.

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

Положим в формуле обращения (4) p = с + гт, тогда exp(pt) = exp(ct) exp(irt). При фиксированном t первый сомножитель постоянен, а второй пробегает единичную окружность на комплексной плоскости бесконечное число раз. С ростом t первый сомножитель и скорость пробегания окружности вторым сомножителем неограниченно возрастают, так что попытка приблизить интеграл в (4) римановыми суммами вряд ли приведет к цели. Например, в простейшем случае для f (t) = 1 имеем F(p) = 1/p, так что при любом с > 0 сомножитель exp(ct) быстро растет с увеличением t, однако, оригинал постоянен и не зависит от с.

Изображение Лапласа функции вида (3) имеет точку ветвления при p = 0. Для устранения многозначности достаточно выбрать одну из ветвей, что делается стандартным образом: приведенные выше преобразования Лапласа не имеют особенностей на комплексной плоскости C \ R- с разрезом вдоль полупрямой

R- = {p £ C : Im(p) = 0, Re(p) < 0}.

С целью уменьшения осцилляций сомножителя exp(iTt) и ограничения скорости роста сомножителя exp(ct) заменим линию интегрирования в (4) контуром L, предложенным в работе [5]:

L=\p\p=z---------2^fL7)--т, У € [-1,1]| •

[ 1 — exp(-2пуг) )

Этот контур состоит из двух симметричных относительно вещественной оси плоскости p ветвей, исходящих из точки p(0) = 1 налево и при у ^ ±1 стремящихся к асимптотам p = ±пг. Для изображения выполнены условия леммы Жордана [6] и кривая L содержит внутри себя особые точки изображения, так что мы имеем

f (t) = J exp(p(y)t)F(p(y))p'y (y) dy

(5)

Для приближенного вычисления этого интеграла применим формулу трапеций.

Заметим, что такой способ страдает дефектом, отмеченным выше — с ростом t при малых y сомножитель exp(p(y)t) неограниченно растет. Поэтому описанный прием пригоден лишь для сравнительно небольших значений t.

Однако это затруднение легко преодолевается предварительной заменой переменной в интеграле (4) вида p = \p\. Так как у изображений нет особых точек справа от мнимой оси, с может быть любым положительным числом. Положим Л = 1/t, тогда экспоненциальный сомножитель под знаком интеграла в (4) становится ограниченным и при y ^ ±1 быстро стремится к нулю.

Подынтегральная функция в представлении

в силу свойств изображения не имеет особенностей как на линии Ь, так и в некоторой полосе, содержащей в себе контур интегрирования Ь. Следовательно, формула трапеций для вычислений интеграла (6) будет давать хорошую точность [7]. Скорость сходимости зависит от ширины полосы и шага численного интегрирования [7].

Пусть теперь контур интегрирования Ь в (6) состоит из нижнего и верхнего берегов разрезов, соединенных окружностью сколь угодно малого радиуса с центром в точке р = 0.

Изображения Г(р) вида (3) фактически зависят от ра, так что положим Г(р) = Ф1_(ра) и введем в рассмотрение функции

Очевидно, F + (t) = F (t) в силу вещественности функции-оригинала. Воспользуемся полученным в работе [8] следующим результатом.

Лемма. Пусть выполнены условия:

(A) F(p) = o(1) при |p| ^ ж, F(p) = o(|p|-1) при |p| ^ 0 равномерно в любом секторе | argp| < п — п, п > п > 0;

(Б) существует е > 0 такое, что для любого ф, удовлетворяющего неравенству п — е < ф < п, справедливы соотношения

(6)

F±(t) = &i(ta exp(±ian)), t> 0.

где a(r) не зависит от p и a(r)exp(-Sr) £ Li(R+) для любого S > 0. Тогда

е xtlmF (t) dt.

1

(7)

При 0 < а +1 < 1, 0 < ^ < 1, в > 0, А > 0 функции (3) удовлетворяют условиям леммы.

В частном случае дробной экспоненты Г(р) = 1/(ра + 1), тогда

F (t) = Im

1

ta sin na

ta exp(—ina) + 1 1 + 2ta eos na, + t2a ’

и формула (7) дает

sin па _xt ta dt

Эа(х) = ------------

п Jо 1 + 2ta cos па + t2a

„ 1 sin па zae-z dz . .

= ж --- ^8

п Jо z2a + 2zaxa cos па + x2a

При x ^ 0 последний интеграл в представлении (8) стремится к величине Г(1 — а), и с учетом формулы Г(а)Г(1 — а) = п/ sin па из представления (8) при x ^ 0 получаем Эа(х) « x^1^^), что совпадает с первым членом ряда (2).

Для приближенного вычисления последнего интеграла можно применить квадратурные формулы типа Гаусса с весом Лагерра zae-z [9]:

оОО ,ь

/ zae-z f (z) dz -Y] Ak f (zk).

J o

к 1

к=1

Эти формулы точны для функций /^) = zm, т = 0,1,..., 2п —1. Однако при убывании а точность формул будет уменьшаться (так как /^) зависит от za).

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

Поэтому для приближенного вычисления интегралов (8) построим специальные обобщенные квадратурные формулы вида

/• ОО п

/ zвв-2/^) dz «У' Ак/^к), в> -1, (9)

к=1

точные для функций /(z) = zam, т = 0, 1 . .., 2п — 1.

Теорема 1. Для того чтобы формула (9) была точна для функций /^) = zam, т = 0, 1 . .., 2п — 1, необходимо и достаточно выполнение двух условий:

1) формула (9) интерполяционная;

2) построенный по узлам формулы (9) многочлен

Áz) — Ц (Z - Za) (10)

k=1

удовлетворяет условиям

О

/ ze e-z LVn(za)zam dz — 0, m — 0,1,...,n — 1. (11)

o

Доказательство. Покажем, что многочлен (10), удовлетворяющий условиям (11), существует и определяется однозначно.

После замены переменной = х условия (11) принимают вид

Xх

(e+1)¡a 1 exp(-x1/a) ^n(x) xm dx — 0, m — 0, 1,...,n — 1.

oo

o

Функция 'ш(х) = х(в+1)/а—1 ехр(—х1/а) обладает свойствами веса на полуоси (0, то), поскольку (в + 1)/а > 0, следовательно, искомый многочлен существует и единствен [9], а его корни, т. е. х^, к = 1, 2,... ,п, попарно различны и положительны. Все коэффициенты формулы положительны. Итак, квадратурная формула типа Гаусса вида (9) существует.

Опишем способ вычисления ее узлов и коэффициентов.

Будем искать многочлен (10) в виде

^и(х) = хп + Ьіхп 1 + • • • + Ъп

с неизвестными коэффициентами Ьк.

Условия (11) приводят к системе линейных алгебраических уравнений

п

Г(в + (к + п — ])а + 1) Ь^ = —Г(в + (к + п)а +1), к = 0,1, ... ,п — 1.

з=1

Ее решение существует и единственно, как показано выше. Далее, находим корни уравнения шп(х) = 0, т. е. числа х^, к = 1, 2,... ,п. Коэффициенты формулы (9) определяем из системы уравнений

п

]Г А к (х аку-1 = Г(в + (з — 1) а + 1), з = 1, 2,... ,п. к =1

Заметим, что узлы и коэффициенты формулы (9) вещественны.

Замечание 1. Описанная схема применима и для других частных случаев ядер (3).

3. Построение комплексных обобщенных квадратурных формул. Для обращения преобразований Лапласа, соответствующих медленно протекающим процессам, целесообразно вместо квадратурных формул наивысшей степени точности (КФНСТ) [4] использовать обобщенные квадратурные формулы наивысшей степени точности (ОКФНСТ) вида

і рС+ІЖ п

— ерр~8ср(р) ¿р ж '^АкФІРк), в>0, (12)

2ПІ О С — ІЖ к =1

точные для <р(р) = р—а,3 = 0,1,..., 2п — 1, или для оригиналов вида ts—1Q2n—l(ta).

Такие формулы были введены в работе [10] и подробно исследованы в статье [11]. В частности, в статье [11] доказана

Теорема 2. Для того чтобы формула (12) была точна для функций ф(р) = р—а° ,з =0, 1, .. ., 2п — 1, необходимо и достаточно выполнение двух условий:

1) формула (12) интерполяционная;

2) построенный по узлам формулы (12) многочлен

п

^п(х) = П (Х — Р—^

=1

удовлетворяет условиям

РС+ІЖ

/ врр—1>шп(р—а)р—ат ¿р = 0, т = 0, 1, ... ,п — 1.

Л С — ІЖ

Доказано [11], что такой многочлен существует и определяется однозначно, а его корни, т. е. узлы ОКФНСТ, удовлетворяют неравенствам Ие (р а) > 0, к = 1, 2,... ,п.

Отметим, что узлы и коэффициенты формулы (12) суть комплексные числа.

Эти формулы оказались весьма эффективными при решении задач линейной вязкоупругости [2].

Вопросы теоретической оценки погрешности КФНСТ и ОКФНСТ рассматривались в работах [12, 13].

Замечание 2. Устойчивость формул (12) по отношению к ошибкам вычисления функции ф(р) определяется суммой модулей коэффициентов А к, которая быстро возрастает при увеличении числа узлов п.

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

Теорема 3. Пусть

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

Ишп^то Ап = ж и |р — ро| < 1о, здесь е„ и Аи зависят от ро. Тогда функция /(£) разлагается в асимптотический ряд

где ^2Р0 означает суммирование по всем особым точкам ро.

Примечание. Если Аи = п — целое неотрицательное число, то 1/Г(—Аи) = 0.

Применим эту теорему к изображениям дробно-экспоненциальной функции и интеграла от нее:

ж < Ао < Ах < ... < Хп < ...,

где а = 1 + а. Следовательно,

п=о

ta+1 J x _n Г(1 — a — an) ’

0 n=0

где x = ta.

5. Обращение преобразования Лапласа с помощью метода Виддера. В

основе этого метода лежит следующая

Теорема 4 ([15]). Если f (t) имеет ограниченную вариацию на отрезке [0, T] при

О

любом T > 0 и интеграл J exp(— pt)f (t) dt = F(p) сходится абсолютно при некотором

0

p, то

Иш tl)lxn+iF(n){Xn) = hf{t _ 0) + /(t + 0)]; n——o n! 2

где xn = (n + 6n)/t, 0 ^ 6n ^ 1, t > 0.

Введем оператор Виддера, полагая вп = 0:

(=). (16)

Очевидно, они существуют для всех п, начиная с некоторого натурального по.

Если функция ](х) непрерывна в точке х = Ь, то при п ^ <х приближение №п(/,Ь) сходится к ](Ь). При этом скорость сходимости будет зависеть от степени гладкости функции и в случае непрерывности второй производной она будет величиной 0(п-1), но дальнейшее увеличение гладкости оригинала не увеличивает скорость сходимости, т. е. метод Виддера быстро насыщаем.

Замечание 3. Функции-оригиналы, соответствующие изображению (3), не ограничены в точке Ь = 0, и непосредственно для них теорема не применима. Можно применить аддитивный метод выделения особенности оригинала в нуле. Поясним его на случае дробной экспоненты (2): введем новую функцию-оригинал

ko 1 (__t l+a ) k

hit) = Эa{t)-ta ]T \ vi-u^V k=0 r((1 + «)(1 + k))

к=0

где ко —наименьшее натуральное число такое, что (1 + а)ко > 1. Эта функция ограничена на всей полуоси и имеет ограниченные производные на любом внутреннем отрезке полуоси, и ее изображение, как легко проверить, имеет вид

*1 (Р) 1

рако (ра _|_ 1)

Если число ако — целое, то эта операция равносильна кратному интегрированию исходного оригинала.

Аналогичное устранение особенности можно осуществить и в других случаях ядер (3), если воспользоваться представлениями оригиналов в виде рядов, приведенных в работе [2].

Существует специальный прием ускорения сходимости метода Виддера.

Выберем к различных натуральных чисел п\ < п2 < ... < пк и положим ]• = ?'/п!,з = 1, 2,... ,к. Составим линейную комбинацию

к

W (n,k,f,t)=Yï, Cjk Wrij (f,t) (17)

j=1

с коэффициентами cjk, определяемыми из системы линейных алгебраических уравнений

к к к к

1 ''ÿ^,cjkdj = 0, cjkdj =0, .. . ^^cjkd- ( ) = 0.

j=1 j=i j=i j=i

Решение этой системы легко находится:

к d

11 у ; j = l,2,...,k. (18)

i=1 j i i=j

Справедлива следующая

Теорема 5 ([16]). Пусть [a1, b1] С [a, b] (0 < a < a1 <b1 <b < ж), и f G C(0, ж), f G С2к([a, b]), тогда при n ^ ж равномерно относительно t G [a1,b1] справедливо равенство

2к

^(W (n, к, f, t) — f (t)) =53 Mmtmf (m)(t) + o(1),

т=к

где Mm — константы, не зависящие от t и f.

Итак, при n ^ ж имеет место равенство W (n,k,f,t) — f (t) = 0(п-к ).

Вычисление Wn(f,t) непосредственно по формуле (16), содержащей производные изображения высокого порядка, затруднительно. Поэтому для преодоления этой трудности был использован предложенный в работе [17] способ, для нашего случая состоящий в следующем: пусть n, m — натуральные числа (m > n) и r G (0,1). Положим

1 ^ (rem (j))-n f n

wn,m(f,t) = -J2 1 (7(1 -rem(j))) , (19)

m 1 - rem(j) \t )

где em(x) = exp(27Г~), <p(p) = pF(p). Предположим, что значения функции <р(р) вычислены с погрешностью, по модулю не превосходящей е, и искомый оригинал ограничен:

|f (t)| ^ M, тогда справедливо неравенство

|Wn,m(f,t) - Wn(f,t)| < д,

где

Mrm е , ч

Д=----------h т-------s—• (20)

1 r'm (1 r)rn

При любом фиксированном n при е ^ 0, т. е. с увеличением точности вычислений функции tp(p), можно выбрать параметры r, m так, что Д ^ 0. Отметим, что значение m следует выбирать из условия, что первое слагаемое справа в оценке (20) приближенно равно второму слагаемому, поскольку увеличение m суммарную погрешность

n

не уменьшает, а лишь увеличивает объем вычислений по формуле (19). Таким образом, для вычислений по формуле (17) вместо Шп^ (/,4) используем приближения к ним ^пз-,„ч (/,4), вычисленные по формуле (19) при подходящем тз.

Рассмотрим теперь следующую задачу: при фиксированном числе слагаемых к в формуле (17) выбрать такие номера приближений по Виддеру из заданного диапазона

к

от П1 до Пк, чтобы коэффициент В = ^2 с^кимел наименьшее возможное значение

з=1

по абсолютной величине (этот коэффициент определяет главный член погрешности приближения (17) к искомому оригиналу). В работе [17] показано, что оно достигается, если использовать приближения Шп с номерами Пк,Пк — 1,... ,Пк — к +1.

Однако при таком выборе чисел ¿1,..., ¿к мы получим наибольшие по модулю значения чисел с3к, и в правой части (17) будут складываться большие числа разных знаков, что может привести к потере точности в вычислениях.

Отсюда возникает следующая задача: для фиксированных значений к,П1,Пк вы-

к

брать числа ¿1,... ,<1к так, чтобы величина У~] |с3к | была минимальна (тогда формула

з=1

(17) наиболее устойчива по отношению к ошибкам в используемых приближениях по Виддеру).

Рассмотрим на отрезке [1/пк, 1/П1] смещенный многочлен Чебышёва

Тк-і(х) = Тк-і

2х--------------------—

П1 Пк

П1

_±_

Пк

к> 1.

Возьмем точки

Х3+1 =

\«1 Пк

008

1 1

н--------1—

П1 Пк

в которых модуль многочлена Тк-і(х) достигает максимума, т. е. |Тк-і(х3-)| = 1, і = 1, 2,..., к. Положим тз = 1/х3', і = 1, 2,..., к, а затем по числам ¿3 = тз/ті, і = 1,2,... ,к, построим числа вида (18). Тогда, как показано в работе [17], величина к „

\С]к\ принимает минимальное значение, равное ( — 1)й_1Та_і(0).

3=1

Числа тз не целые, и поэтому следует в качестве искомых номеров приближений по Виддеру взять ближайшие целые к ним, т. е. положить П3 = [тз + 1/2], і = 1, 2, ... ,к. Такие номера будем называть чебышёвскими.

Замечание 4. Все описанные выше методы обращения преобразования Лапласа были реализованы в виде программ на различных языках программирования. Полученные результаты численных экспериментов были сопоставлены с табличными данными в цитированных работах, что позволяет сделать вывод о высокой эффективности рассмотренных алгоритмов.

1

2

Литература

1. Работное Ю. Н. Элементы наследственной механики твердых тел. М., 1977. 384 с.

2. Екельчик В. С., Рябое В. М. Об использовании одного класса наследственных ядер в линейных уравнениях вязкоупругости // Механика композитных материалов. 1981. №3. С. 393404.

3. Работнов Ю. Н., Паперник Л. Х., Звонов Е. Н. Таблицы дробно-экспоненциальной функции отрицательных параметров и интеграла от нее. М., 1969. 132 с.

4. Крылов В. И., Скобля Н. С. Методы приближенного преобразования Фурье и обращения Лапласа. М., 1974. 224 с.

5. Talbot A. The accurate numerical inversion of Laplace transform // J. Inst. Maths. Applies. 1979. Vol. 23. P. 97-120.

6. Лаврентьев М.А., Шабат Б. В. Методы теории функций комплексного переменного. М., 2002. 688 с.

7. Самокиш Б. А. Замечание о вычислении определенных интегралов // Методы вычислений. Вып. 2. Л., 1963. С. 45-49.

8. Bobylev A. V., Cercignani C. The inverse Laplace transform of some analytic functions with an application to the eternal solutions of the Boltzmann equation // Applied Mathematics Letters. 2002. Vol. 15. P. 807-813.

9. Мысовских И. П. Лекции по методам вычислений. СПб., 1998. 472 с.

10. Рябов В. М. О многочленах, возникающих при численном обращении преобразования Лапласа // Методы вычислений. Вып. 12. Л., 1981. С. 46-53.

11. Рябов В. М. Свойства квадратурных формул наивысшей степени точности, применяемых для обращения преобразования Лапласа // Журнал вычисл. матем. и матем. физ. 1989. Т. 29. №7. С. 1083-1087.

12. Матвеева Т. А., Рябов В. М. Об оценках погрешности квадратурных формул численного обращения преобразования Лапласа // Вестн. С.-Петерб. ун-та. Сер. 1. 2000. Вып. 4. С. 7-11.

13. Матвеева Т. А., Рябов В. М. Обоббщенные квадратурные формулы численного обращения преобразования Лапласа // Вестн. С.-Петерб. ун-та. Сер. 1. 2002. Вып. 4. С. 27-33.

14. Диткин В. А., Прудников А. П. Операционное исчисление. М., 1975. 407 с.

15. Widder D. V. The Laplace transform. Oxford, 1946. 406 p.

16. May C. P. Saturation and inverse theorems for combinations of a class of exponential-type operators // Canad. J. Math. 1976. Vol. 29. N6. P.1224-1250.

17. Рябов В. М. О некоторых задачах, возникающих при обращении преобразования Лапласа // Методы вычислений. Вып. 16. Л., 1991. С. 59-67.

Статья поступила в редакцию 17 февраля 2009 г.