Диссертация (1150465), страница 9
Текст из файла (страница 9)
Скорость его сходимости можно увеличить,построив линейные комбинации fn (t), вычисленные для различных n [37]. Однако для этого могут потребоваться формулы (56) с большим числом узлов m.Рассмотрим алгоритм для их построения.Узлы формулы (56) совпадают с корнями многочленов Лагерра, т. е.Lm (pk ) = 0,k = 1, m.(57)Будем находить корни уравнения (57) методом Ньютона:pj+1= pjk − Lm (pjk )/L0m (pjk ),kj = 0, 1, .
. .(58)Начальные приближения ко всем корням многочленов (55), для которых метод(58) сходится, берем из приближенных равенств [60] (нам нужен только случайα = 0, т. к. Lm (p) = Lm (p ; 0))(1 + α)(3 + 0.92α),1 + 2.4m + 1.8α15 + 6.25αp2 − p1 ≈,1 + 0.9α +2.5m11.26kαpk+2 − pk+11 + 2.55k≈+,pk+1 − pk1 + 0.3α1.9k1 + 3.5kp1 ≈k = 1, m − 2.После нахождения с требуемой точностью всех узлов вычисляем коэффициентыформулы (56):Ak =pk,(L0m (pk ))2k = 1, m.71Для вычисления многочленов Лагерра следует использовать рекуррентное соотношение(n + 1)Ln+1 (x ; α) = (α + 2n + 1 − x)Ln (x ; α) − (α + n)Ln−1 (x ; α),n = 1, 2, .
. . ,L0 (x ; α) = 1,L1 (x ; α) = 1 + α − x,а значения их производных исключать с помощью соотношенияxL0n (x ; α) = nLn (x ; α) − (α + n)Ln−1 (x ; α).Так, получим выражение для L0m (pk )L0m (pk ) =−mLm−1 (pk ).pkВ [44] был предложен другой подход к построению функции–оригинала,использующий значения изображения в равноотстоящих точках вещественнойполуоси p > 0. Исследование этого метода было проведено в работах [35], [36],там же приведены числовые параметры конкретных частных случаев общейсхемы. Многочлены Лагерра, которые используются с целью ускорения сходимости метода, изучались в работах [2], [52].723.3.Аддитивный метод вычисления дробно–экспоненциальнойфункцииКак отмечалось ранее изображения по Лапласу функции Эα (t) и интеграRtла от нее 0 Эα (τ )dτ равны, соответственно,1,pa + 11,p(pa + 1)a = 1 + α.Для вычисления функции F1 (α, x) воспользуемся аддитивным методом.Рассмотрим его сущность сначала на примере изображения самой дробно–экспоненциальной функции —1pa +1 .Заметим, что изображение представимо в видеряда∞X (−1)k1=.pa + 1pa(k+1)k=0Такое представление следует из следующих выкладок, где мы воспользовалисьформулой суммы для бесконечно убывающей геометрической прогрессии:k∞ 1111 X −1=·.=1papapa 1 + papa + 1k=0Разобьём исходную сумму на двеkX∞∞0 −1XX(−1)k(−1)k(−1)k=+,a(k+1)a(k+1)a(k+1)pppk=0k=0k=k(59)0где k0 ∈ N определяется из условия, которое будет выведено позже.Для вычисления второй суммы в правой части равенства (59) воспользуемся формулой суммы бесконечно убывающей геометрической прогрессии сознаменателем q = − p1a .
В результате получим следующее выражение:∞X(−1)k1k0.=(−1)a(k+1)(pa + 1)pak0pk=k0Полученное изображение, с другой стороны, есть разность между преобразованием Лапласа дробно–экспоненциальной функции и первой суммы из (59):kX∞0 −1X(−1)k(−1)k= L(Эα (t)) −.a(k+1)a(k+1)ppk=kk=00(60)73Отсюда видно, что функция–оригинал для этого изображения есть разностьмежду дробно–экспоненциальной функцией и функцией–оригиналом для суммы в правой части равенства (60). То есть значение дробно–экспоненциальнойфункции находится как сумма двух оригиналов.
Рассмотрим задачу обращениядля обоих изображений в отдельности.В силу линейности преобразования Лапласа оригиналом для изображенияkX0 −1k=0(−1)kpa(k+1)будет сумма функций–оригиналов для каждого члена ряда, при условии ak0 ≥−α, k0 ∈ N, чего можно добиться соответствующим выбором параметра k0 .Таким образом получимL−1kX0 −1k=0k(−1)pa(k+1)!=kX0 −1k=0(−1)k · tak+a−1,Γ(ak + a)где ak0 > −α, k0 ∈ N.1Теперь рассмотрим задачу обращения изображения (−1)k0 (pa +1)pak0 .
Дляэтого воспользуемся ОКФНСТ, которые были рассмотрены в Главе 1 Раздел1.3. Для нашего случая функция ϕ(p) = ps F (p) в (12) имеет вид(−1)k0pa +1 ,такимобразом параметр s = ak0 . Используя теорему 3 из Главы 1 Раздел 1.3, былинайдены узлы pk и коэффициенты Ak ОКФНСТ, тем самым построен методдля нахождения значений оригинала для изображения(−1)k0.(pa + 1)pak0В результате получим конечную формулу для нахождения значений дробно–экспоненциальной функции:Эα (t) ≈ tak0 −1nXk=1где ϕ(p) =(−1)k0pa +1 ,Ak ϕ(pk /t) +kX0 −1k=0(−1)k · tak+a−1,Γ(ak + a)k0 ∈ N и находится из условия ak0 > −α.74Описанный алгоритм был реализован в среде Maple для нахождения значений функции F1 (α, x).
Также была посчитана оценка погрешности вычисления самой функции F1 (α, x) и оценки погрешности для ОКФНСТ применительно к нашему случаю, которые были рассмотрены в статье [24]|εn (ta )| 6 Maρs · σn (ta , ρ),σn (ta , ρ) = ∞Xm=2nMaρs =1−Γ(s + am)12πnXk=1!2 1/2 ta 2mmAk µk· ,ρ2πZ(61)2(62)1/2|Pas (ρ exp(iθ))| dθ,(63)0aгде µk = p−ak , 0 < t < ρ 6 1.Покажем, как будет вычисляться оценка погрешности ОКФНСТ для нашего случая. Понятно, что σn (ta , ρ) останется такой же, изменится лишь формула для вычисления Maρs . Функция Pas (t) = ϕs (t−1/a ), где s = ak0 , есть(−1)k0 · tPas (t) =.1+tТем самым имеем(−1)k0 · ρ exp(iθ)Pas (ρ exp(iθ)) =.1 + ρ exp(iθ)Подставим выражения для функции Pas (t) в формулу (63) и получим оценкусверху для величины Maρs .Maρs =12πZ2π2|Pas (ρ exp(iθ))| dθ1/26012πρ2·2π (1 − ρ)21/26ρ.1−ρВ результате оценка погрешности ОКФНСТ будет|εn (ta )| 6∞Xρ·1−ρm=2naгде µk = p−ak , 0 < t < ρ 6 1.1−Γ(s + am)nXk=1!2 1/2 ta 2mmA k µk· ,ρ75Оценка погрешности (δ) вычисления функции F1 (α, x) зависит от оценкипогрешности ОКФНСТ и вычисляется по следующей формулеts−1 εn (ta )δ=.tαОписанный метод был реализован в виде программы в математическомпакете Maple.
Текст программы приведён в Приложении Программа 8. Здесьприведем результаты вычисления функции F1 (α, x).Результаты вычисления F1 (α, x) при α = −0.75,x = 0.1.Таблица 21nk0F1F1табл.δ530.226657 0.2266573570.226657 0.2266573 6.673406 · 10−171030.226657 0.2266573 2.764635 · 10−221070.226657 0.2266573 7.942103 · 10−289.80489 · 10−12Результаты вычисления F1 (α, x) при α = −0.75,x = 1.1.Таблица 22nk0F1F1табл.δ530.057488 0.5749 · 10−1 4.063073 · 10−8570.057488 0.5749 · 10−11030.057488 0.5749 · 10−1 3.92322 · 10−181070.057488 0.5749 · 10−1 1.62817 · 10−193.93269 · 10−9Как видно из Таблиц 21, 22, оценка погрешности имеет порядок 10−10 , вто время как табличные значения функции F1 (α, x) из книги [30] даны с точностью до четырёх знаков после запятой. Тем самым значения, получаемые спомощью применения аддитивного метода, являются более точными по сравнению с табличными значениями.763.4.Метод обращения преобразования Лапласа с помощьюразложения оригинала в обобщенные степенные рядыДля больших значений аргументов целесообразно использовать другойподход к вычислению, который основан на следующей теореме (см.
[9]):Теорема 11. Пусть1f (t) = limω→∞ 2πiσ+iωZept F (p) dp ,σ−iωпричем функция F (p) удовлетворяет условиям леммы Жордана, имеет конечное число особых точек и все особые точки являются полюсами или точками ветвления, и в окрестностях особых точек p = p0 с наибольшей вещественной частью функция F (p) разлагается в ряды вида∞Xcν (p0 )(p − p0 )λν ,−∞ < λ0 < λ1 < . . . < λn < . . . ,ν=0limn→∞ λn = ∞ и |p − p0 | < l0 , здесь cν и λν зависят от p0 .
Тогда функция f (t)разлагается в асимптотический рядf (t) ≈Xp0гдеPp0p0 te∞Xcν (p0 ) −λν −1t,Γ(−λ)νν=0означает суммирование по всем особым точкам p0 .П р и м е ч а н и е. Если λν — целое неотрицательное число, то1/Γ(−λν ) = 0.Применим эту теорему к изображениям дробно–экспоненциальной функции и интеграла от нее:∞X1=(−1)n pan ,a1+pn=0∞X1pa−111a−1= − a= −p(−1)n pan ,ap(p + 1) p p + 1 pn=077где a = 1 + α. Следовательно,∞−n−1X1n xF1 (α, x) = α Эα (t) =(−1),tΓ(−an)n=0F2 (α, x) =1Zttα+10∞1 X (−1)n x−n,Эα (τ )dτ = 1 −x n=0 Γ(1 − a − an)где x = ta .Для ядра Гаврильяка–Негами получим следующее разложение:∞X1γk pak =F (p, a, b, c, d) = a=c(p + b) + dk=01bc c− cpa +c2b + d (b + d) bc−b c(c − 1)/2(bc + d)b2 + (bc )2 c2 /(bc + d)2 b2 2a+p − ...,bc + dиз которого по теореме 11 получаем асимптотическое разложение при большихt:N (t, a, b, c, d) =∞Xk=1γk.Γ(−ak)tak+1Для вычисления при больших t функцииZ tN (x, a, b, c, d) dx,0имеющей изображение F (p, a, b, c, d)/p, аналогично получаем формулуZtN (x, a, b, c, d) dx =0∞Xk=0γk.Γ(1 − ak)takОписанный метод был реализован в виде программы в математическомпакете Maple.
Текст программы приведён в Приложении Программа 9. Здесьприведем результаты вычисления функций F1 (α, x) и F2 (α, x), полученные спомощью асимптотических формул.78Результаты вычисления F1 (α, x) и F2 (α, x) при α = −0.65.Таблица 23nxF1F1табл.F2F2табл.10 900 3.117 · 10−7 3.107 · 10−7 1.110 · 10−3 1.109 · 10−310 600 7.009 · 10−7 6.993 · 10−7 1.665 · 10−3 1.665 · 10−310206.024 · 10−4 6.024 · 10−4 4.824 · 10−2 4.824 · 10−21041.222 · 10−2 1.222 · 10−2 2.101 · 10−1 2.101 · 10−11023.809 · 10−2 3.731 · 10−2 3.589 · 10−1 3.590 · 10−12023.755 · 10−2 3.731 · 10−2 3.590 · 10−1 3.590 · 10−1Результаты вычисления F1 (α, x) и F2 (α, x) при α = −0.85.Таблица 24nxF1F1табл.F2F2табл.10 900 1.661 · 10−7 1.655 · 10−7 1.110 · 10−3 1.110 · 10−310 600 3.735 · 10−7 3.733 · 10−7 1.664 · 10−3 1.664 · 10−310203.099 · 10−4 3.099 · 10−4 4.785 · 10−2 4.785 · 10−21045.684 · 10−3 5.666 · 10−3 2.038 · 10−1 2.038 · 10−12021.627 · 10−2 1.626 · 10−2 3.435 · 10−1 3.435 · 10−1Как видно из Таблиц 23, 24 при малых значениях x < 2 точность плохая.
При увеличении числа слагаемых в бесконечной сумме точность метода неулучшается. Для больших x > 2 для достижения хороших результатов вполнедостаточно число слагаемых равное 10. Причём, чем меньше α, тем лучше точность асимптотических формул.79ЗаключениеПеречислим основные результаты работы:1. Исследованы различные методы обращения преобразования Лапласа ипредложены алгоритмы по применению методов обращения к вычислениюфункций специального вида.2. Рассмотрены задачи линейной вязкоупругости и применение преобразования Лапласа для их решения.3. Исследованы свойства ядер, которые могут быть выбраны в качестве ядерползучести и релаксации в соотношении Больцмана–Вольтерра; изученыих свойства.4.
Исследованы методы обращения преобразования Лапласа в предположении, что заданное изображение F (p) искомой функции–оригинала фактически зависит от 1/pa , где a — произвольное положительное число изинтервала (0, 1); получены новые формулы, обладающие большей точностью по сравнению с известными для определенного класса изображенийи имеющие большое прикладное значение.5.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
