Диссертация (Численное обращение интегрального преобразования Лапласа функций специального вида), страница 4

. , n) — узлы икоэффициенты ОКФНСТ и предположим, что |p−ak | < 1, k = 1, 2, . . . , n. Тогдасоответствующая функция–оригинал представима в виде!np XkAk ϕsf (t) = ts−1+ εn (ta ) ,tk=1и для любых чисел t и ρ, удовлетворяющих неравенству 0 < ta < ρ 6 r,справедлива оценка|εn (ta )| 6 Masρ σn (ta , ρ),гдеMasρ =σn (ta , ρ) = ∞Xm=2n12πZ2π|Pas (ρ exp(iθ))|2 dθ 12,01−Γ(s + am)nX!2Ak µmk1/2|ta /ρ|2m ,k=1µk = p−ak .Заметим, что оценка погрешности, доказанная в теореме 5, состоит из двухсомножителей: Masρ , зависящего лишь от изображения и σn (ta , ρ), зависящегоот параметров ОКФНСТ.Метод применения ОКФНСТ для обращения преобразования Лапласа реализован в виде программы в математическом пакете Maple. Текст программыприведён в Приложении Программа 2.

В Таблицах 3, 4 приведены результатывычисления функции F1 (α, x).25Результаты вычисления F1 (α, x) при α = −0.4.Таблица 3nxF1F1табл.20 0.2 4.9091 · 10−1 4.9091 · 10−120 0.8 2.1674 · 10−1 2.1674 · 10−1201.7110 · 10−1 1.7110 · 10−11Результаты вычисления F1 (α, x) при α = −0.65.Таблица 4nxF1F1табл.2019.1239 · 10−2 9.124 · 10−21041.2219 · 10−2 1.222 · 10−22041.2219 · 10−2 1.222 · 10−220 400 1.5759 · 10−6 1.576 · 10−6Как видно из результатов, ОКФНСТ обладают более высокой точностьювычислений по сравнению с КФНСТ. При этом количество узлов n в ОКФНСТбралось равным 20. Данный метод может быть применён для вычисления значений функции–оригинала как при x > 1, так и при 0 < x < 1.Теперь приведём результаты вычисления функции ползучести (пространственная переменная x опущена)Ztε(t) = 1 + 4Э−0.8 (τ )dτ0по её изображениюε(p) =14+.p p(p0.2 + 1)В Таблицах 5–7 приведены результаты вычислений функции ε(t) с помощьюКФНСТ и ОКФНСТ.

Для её вычисления использовались также Программы 1,2 из Приложения с изменённым выражением (изображения по Лапласу) для26подсчёта квадратурной суммы в конце программы.Результаты вычисления ε(t) при α = −0.8,t = 0.001.Таблица 5nε(t)КФНСТ ε(t)ОКФНСТ ε(t)точное101.869591.866961.86696151.865441.866961.86696201.867991.866961.86696Результаты вычисления ε(t) при α = −0.8,t = 1.Таблица 6nε(t)КФНСТ ε(t)ОКФНСТ ε(t)точное103.122883.115603.11559153.111193.115603.11559203.118673.115603.11559Результаты вычисления ε(t) при α = −0.8,t = 106 .Таблица 7nε(t)КФНСТ ε(t)ОКФНСТ ε(t)точное104.799094.793474.79352154.789314.793474.79352204.796824.793474.79352Как видно из приведённых выше таблиц, ОКФНСТ позволяют находитьзначение функции–оригинала с большей точностью по сравнению с КФНСТ.Причём наилучшие результаты достигаются при вычислениях в точках 0<t<1.Таким образом, метод, основанный на применении ОКФНСТ для вычисленияфункции–оригинала по её изображению по Лапласу, зависящего от 1/pa , где0<a<1, является предпочтительным.

При этом число узлов в ОКФНСТ достаточно брать равным n = 10.27Метод построения ОКФНСТ для обращения интегрального преобразования Лапласа рассматривался и изучался в работах [26], [43].281.4.Обращение преобразования Лапласа с помощью методаВиддераНаряду со специальными методами решения интегро–дифференциальныхуравнений с дробными производными могут быть применены также методыВиддера, деформации контура и другие.

Так, рассмотрим метод обращения, воснове которого лежит следующаяТеорема 6 ([63]). Если f (t) имеет ограниченную вариацию на отрезкеR∞[0, T ] при любом T > 0 и интеграл exp(−pt)f (t) dt = F (p) сходится абсо0лютно при некотором p, то1(−1)n n+1 (n)limxn F (xn ) = [f (t − 0) + f (t + 0)],n→∞n!2где xn = (n + θn )/t, 0 6 θn 6 1.Введём оператор Виддера (Wn ), полагая θn = 0:nn+1 (n) n Wn (f, t) = (−1)F,n! tn+1tnn = 1, 2, . .

.(16)Если функция f (x) непрерывна в точке x = t, то при n → ∞ приближениеWn (f, t) сходится к f (t). Cкорость сходимости будет зависеть от степени гладкости функции и в случае непрерывности второй производной она будет величиной O n1 , но дальнейшее увеличение гладкости оригинала не увеличиваетскорость сходимости, т. е.

метод Виддера быстро насыщаем.Существует специальный метод для ускорения сходимости метода Виддера [42]. Выберем k различных натуральных чисел n1 < n2 < · · · < nk и положимdj = nj /n1 , j = 1, 2, . . . , k. Составим линейную комбинациюW (n, k, f, t) =kXcjk Wnj (f, t)(17)j=1с коэффициентами cjk , определяемыми из системы линейных алгебраических29уравненийkXkXcjk = 1,j=1kXcjk d−1j = 0,j=1cjk d−2j= 0,kX...j=1−(k−1)cjk dj= 0.j=1Решение этой системы легко находится:cjk =kYi=1i6=jdj,dj − dij = 1, 2, . . .

, k.(18)Справедлива следующаяТеорема 7 ([59]). Пусть [a1 , b1 ] ⊂ [a, b] (0 < a < a1 < b1 < b < ∞), и f ∈C(0, ∞), f ∈ C 2k ([a, b]), тогда при n → ∞ равномерно относительно t ∈[a1 , b1 ] справедливо равенствоkn (W (n, k, f, t) − f (t)) =2kXMm tm f (m) (t) + o(1),m=kгде Mm — константы, не зависящие от t и f .Итак, при n → ∞ имеет место равенство W (n, k, f, t) − f (t) = O(n−k ).Вычисление Wn (f, t) непосредственно по формуле (16), содержащей производные изображения высокого порядка, затруднительно. Поэтому для преодоленияэтой трудности был использован путь, предложенный в работе [32] и для нашего случая состоит в следующем: пусть n, m — натуральные числа (m > n) иr ∈ (0, 1). Положимm1 X (rem (j))−n nWn,m (f, t) =ϕ(1 − rem (j)) ,m j=1 1 − rem (j)t(19)где em (x) = exp(2π ixm ), ϕ(p) = pF (p).Покажем, что значение Wn (f, t) можно заменить с необходимой точностьюзначением величины Wn,m (f, t).

Поскольку в работе [32] опущен вывод формулы30(19), поэтому приведём его здесь. Рассмотрим функцию11−zG(z, t) = F, t > 0.ttПри фиксированном t > 0 функция G(z, t) регулярна в круге |z| < 1 и разлагается в сходящийся в этом круге рядG(z, t) =∞Xan (t)z n .n=0Коэффициенты an (t) можно вычислять, применяя формулу КошиZ1an (t) =2πi1G(z, t)dz=z n+12πiG(reiϕ , t)ireiϕdϕ(reiϕ )n+10|z|=r1=2πZ2πZ2π(reiϕ )−n G(reiϕ , t)dϕ,0где r — любое положительное вещественное число, r < 1.Для вычисления последнего получившегося интеграла применим формулуправых прямоугольников с шагом h = 2π/m, ϕj = jh, где j = 1, 2, .

. . , m.Получимm1 X iϕj −nan (t) ≈(re ) G(reiϕj , t).m j=1Как показано в работе [42], tWn (f, t) = an.n mt1 X iϕj −n(re ) G(reiϕj , t/n) = Wn,m (f, t),Wn (f, t) = an≈nm j=1таким образом получим формулу (19).Предположим, что значения функции ϕ(p) вычислены с погрешностью,по модулю не превосходящей ε, и искомый оригинал ограничен: |f (t)| 6 M ,тогда справедливо неравенство|Wn,m (f, t) − Wn (f, t)| 6 ∆,31гдеM rmε∆=+.(20)m1−r(1 − r)rnПри любом фиксированном n при ε → 0, т. е. с увеличением точности вычислений функции ϕ(p), можно выбрать параметры r, m так, что ∆ → 0. Отметим,что значение m следует выбирать из условия, что первое слагаемое справа воценке (20) приближенно равно второму слагаемому, поскольку увеличение mсуммарную погрешность не уменьшает, а лишь увеличивает объем вычисленийпо формуле (19).Следовательно возникает задача нахождения оптимальных значений параметров r, m, при которых ∆ → 0.

Для их нахождения запишем условие (20)в видеε≈ δ,(1 − r)rnrmM≈ δ.1 − rmЗаданными и известными величинами в условиях, приведённых выше, являются ε — точность вычислений (либо машинная точность, либо точность вычисления функции ϕ(p)), M — величина, ограничивающая модуль функции–оригинала. Задаём параметры n и δ (необходимая точность вычисления пометоду Виддера) и из первого условия находим значение r. Затем найденноезначение r подставляем во второе условие и находим величину m.

Используянайденные значения величин r, m, находились коэффициенты cjk по описаннымвыше формулам.Таким образом для вычислений по формуле (17) вместо Wnj (f, t) используем приближения к ним Wnj ,mj (f, t), вычисленные по формуле (19) при подходящем mj .Рассмотрим теперь следующую задачу: при фиксированном числе слагаемых k в формуле (17) выбрать такие номера приближений по Виддеру изkPзаданного диапазона от n1 до nk , чтобы коэффициент B =cjk d−kj имел наиj=1меньшее возможное значение по абсолютной величине (этот коэффициент опре-32деляет главный член погрешности приближения (17) к искомому оригиналу).В работе [32] показано, что оно достигается, если использовать приближенияWn с номерами nk , nk − 1, .

. . , nk − k + 1.Однако при таком выборе чисел d1 , . . . , dk мы получим наибольшие по модулю значения чисел cjk , и в правой части (17) складываются большие числаразных знаков, что может привести к потере точности в вычислениях. Отсюдавозникает следующая задача: для фиксированных значений k, n1 , nk выбратьkPчисла d1 , . . . , dk так, чтобы величина|cjk | была минимальна (тогда формулаj=1(17) наиболее устойчива по отношению к ошибкам в используемых приближениях по Виддеру).Рассмотрим на отрезке [ n1k , n11 ] смещенный многочлен Чебышева!2x − n11 − n1kTek−1 (x) = Tk−1, k > 1.11−n1nkВозьмём точкиxj+1вкоторыхмодуль|Tek−1 (xj )| = 1,числам d¯j =1=211−n1 nkмногочленаjπcosk−1Tek−1 (x)j = 1, 2, .