CC BY
12_ВЕСТНИК ПЕРМСКОГО УНИВЕРСИТЕТА_
2024 • Математика. Механика. Информатика« Вып.1(64)
«Математика»
Научная статья УДК 519.6
DOI: 10.17072/1993-0550-2024-1-5-14
Модифицированная формула Герасимова-Капуто
Н.К. Волосова1, К.А. Волосов2, А.К. Волосова2, М.И. Карлов3, Д.Ф. Пастухов4, Ю.Ф. Пастухов4
Московский государственный технический университет (МГТУ) им. Н.Э. Баумана, Москва, Россия 2 Российский университет транспорта (МИИТ), Москва, Россия 3Московский физико-технический университет (МФТИ), Москва, Россия 4Полоцкий государственный университет, Новополоцк, Республика Беларусь
Автор, ответственный за переписку: Дмитрий Феликсович Пастухов, [email protected]
Аннотация. В работе впервые получены модифицированные формулы Герасимова-Капуто. Модифицированные формулы учитывают значение производной функции в нуле с порядком на единицу меньше, чем порядок производной, стоящей под знаком интеграла Герасимова-Капуто. Без учета нового слагаемого в формулах Герасимова-Капуто не всегда корректно вычисление дробной производной на интервале любого порядка и для любой функции. В работе также описан простой численный алгоритм с квадратурной формулой Гаусса, позволяющей вычислять дробную производную с двойной точностью. Составлены таблицы дробной производной для функций синуса и косинуса. Причем первая половина таблиц (в интервале порядка (0,1)) и вторая половина таблиц (в интервале порядка (1,2)) получена программами по разным алгоритмам. В работе достигнута абсолютная погрешность вычисления дробной производной 10-15.
Ключевые слова: численные методы; дробная производная Герасимова-Капуто; дробная производная Капуто
Для цитирования: Волосова Н.К., Волосов К.А., Волосова А.К., Карлов М.И., Пастухов Д.Ф., Пастухов Ю. Ф. Модифицированная формула Герасимова-Капуто // Вестник Пермского университета. Математика. Механика. Информатика. 2024. Вып. 1(64). С. 5-14. DOI: 10.17072/ 1993-0550-2024-1-5-14.
Статья поступила в редакцию 19.01.2024; одобрена после рецензирования 14.02.2024; принята к публикации 15.03.2024.
«Mathematics»
Research article
Modified Gerasimov-Caputo Formula
N.K. Volosova1, K.A. Volosov2, A.K. Volosova2, M.I. Karlov3, D.F. Pastuhov4, Yu.F. Pastuhov4
:Bauman Moscow State Technical University (BMSTU), Moscow, Russia 2Russian University of Transport (RUT MIIT), Moscow, Russia 3Moscow University of Physics and Technology (MIPT), Moscow, Russia 4Polotsk State University, Novopolotsk, Republic of Belarus
Corresponding author: Dmitriy F. Pastukhov, [email protected]
Данная работа О 2024 Волосова Н.К., Волосов К.А., Волосова А.К., Карлов М.И., Пастухов Д.Ф., Пастухов Ю.Ф. распространяется под лицензией СС BY 4.0. Чтобы просмотреть копию этой лицензии, посетите http://creativecommons.org/licenses/by/4.0/
Abstract. In this work, modified Gerasimov-Caputo formulas were obtained for the first time. The modified formulas take into account the value of the derivative of the function at zero with an order of one less than the order of the derivative under the sign of the Gerasimov-Caputo integral. Without taking into account the new term in the Gerasimov-Caputo formulas, it is not always possible to calculate the fractional derivative on any order interval and for any function. The paper also describes a simple numerical algorithm with the Gaussian quadrature formula, which allows one to calculate the fractional derivative with double precision. tables of fractional derivatives for the sine and cosine functions have been compiled. Moreover, the first half of the tables (in the interval of order (0,1)) and the second half of the tables (in the interval of order (1,2)) were obtained by programs using different algorithms. In the work, an absolute error in calculating the fractional derivative of 10-15 was achieved.
Keywords: numerical methods; Gerasimov-Caputo fractional derivative; Caputo fractional derivative
For citation: Volosova, N.K., Volosov, K.A., Volosova, A.K., Karlov, M.I., Pastukhov, D.F., Pastukhov, Yu.F.
(2024), "Modified Gerasimov-Caputo Formula", Bulletin of Perm University. Mathematics. Mechanics. Computer
Science, no. 1(64), pp. 5-14. DOI: 10.17072/1993-0550-2024-1-5-14.
The article was submitted 19.01.2024; approved after reviewing 14.02.2024; accepted for publication 15.03.2024.
Введение
В работе [1] показано, что в задачах аномальной диффузии необходимо использовать производные дробного порядка, принимающего значения на интервале (0,2). Другой пример, в задачах механики поток газа Трикоми на звуковой лини прямо пропорционален производной порядка 2/3 от функции тока [2]. Как известно, производная Гераси-мова-Капуто определяется [1], [2], [3], [4] интегральной формулой (1) .
В работе показано, что для явного вычисления производной Герасимова-Капуто необходимо еще одно слагаемое, которое зависит от значения функции и(0), переменной t и порядка дробной производной а на интервале (0,1). Данное слагаемое необходимо для согласования производной функции целого порядка и предельного значения производной дробного порядка, когда ее порядок становится целым. Для некоторых функций новое дополнительное слагаемое в точке 1=0 не требуется, например, для функции и(1)=8т(1) в интервале а = (0,1). В работе максимально упрощен алгоритм вычисления дробных производных квадратурной формулой Гаусса. В программу вводится исходная функция, первая или вторая ее производные. Гладкости первого или второго порядка используемой функции требует также формула-определение (1).
В работах [5], [6], [10] рассмотрены разностные уравнения и аналог задачи оптимального управления Л.С. Понтрягина с производными дробного порядка.
Постановка задачи
Пусть заданная функция и(х,1) является достаточно гладкой по переменной t, и(х, ?) е Сп (0, ?), тогда производная функции целого порядка п по переменной t является непрерывной и интеграл в формуле (1) сходится, так как в особой точке г ^ / знаменатель дроби пропорционален 1 /(/ -г)а~"+1 учитывая -1 <а -п +1 <0 .
Определение 1. Производной дробного порядка а>0 Герасимова-Капуто от функции двух переменных u(x,t) (по переменной 0 [1], [2], [3] называется функция
(Diu )(x, t )=
1 г д"u(x,z)
Г(n - а)
dz
dzn
(t -z)a
(1)
i = [а] +1,n e N,а e R,а - n +1 = {а}е (0,l).
Для простоты рассмотрим частный случай дробной производной для функции одной переменной u(t), так как в нашей задаче переменная х не используется. Для двух интересующих нас интервалов порядка производной с учетом формулы (1) получим:
1). 0 <а< 1, n = 1
(Dku )(t )=
1
г du(z) dz , „ , ,.„4 J—-—, n = 1,0 < а < 1 (2)
Г(1 - а )J dz (t -z)
2). 1 < а < 2,n = 2 а V\ 1 fd 2u(z) dz „, „
)(t '=TfOiJ ^ ырn=21 < а <2 .(3)
0
Гамма-функция в формулах (2), (3) вычисляется по формуле Эйлера (4):
WUk)
1
'rdu(r) dr
Г(1 -a)|J dr (t -r)c
- + u(0)gi(t,a)
D+A )=^ [J^TAt - du(0)g2M)
Г(2-a)IJ dr2 (t-rf
dt
для a e (0,1):
u(0)
1 'rdu(r) dr Г(1 -a)0 dr (t-r)a
Г(a)= J ta-1etdt, t e [0, <»),a> 0 . (4) = Г(1 -a)
+ lim
ta s^+0
( и Л t-s
u(t - s)
w
u(r)
- at , u(r) , dr
0 (t -r)
.(9)
Поскольку дробные производные в формулах (2), (3) определены с точностью до множителя Г (i -а), Г (2 -а), то сложность вычисления производной дробного порядка заключается в сложности численного алгоритма для соответствующего интеграла. В языке FORTRAN вызовом функции dgamma(a) Г (а) вычисляется с относительной точностью 10-i5.
Очевидно, что производная дробного порядка а в (2), (3) должна переходить в производную целого порядка n-1 в случае а ^ n -1. Если в формуле (2) а ^ n -1 = 0 , то
В формуле (9) под знаком предельного перехода два слагаемые, зависящие от переменной е, стремятся к бесконечности, и именно их разность может дать конечное число. Уточняя формулу (7), получим
1 ('rdu(r) dr u(0)^
Г(1 -a)[° dr (t -r)a + Ia
Й+Л) = ^ I J.
V 0
Преобразуем также (3):
ae (0,1)(. (10)
1 t d2u(x) dr
J T fn-=
te"u)(') = 4\~¡~77~TY=0 = urn -u(0). (5) Г(2""i
u (0)
—-Г + lim
ta-1 s-^+0
Г(1)J0 dr (t-r)
Аналогично, если в формуле (3) а ^ n -1 = 1, то
1 fd2u(r) dr d , „ d
(01 )=^ Í -г1^=- dtu(0). (6)
Рассмотрим пример:
п (5)
u(t) = sin(t) : u(t) = sin(t) - sin(0) = sin(t)
, (6)
D%+u(t) = sin (t) - sin (0) = cos(t) -1 Ф cos(t).
Получаем противоречие для производной синуса первого порядка и формулой (6) на интервале а е (1,2), хотя противоречия для производной синуса нуль порядка на интервале а е (0,1) и формулой (5) нет.
Чтобы устранить противоречие, в формулу (2) нужно добавить слагаемое вида
(u(tW
(s)a-1
- (a-1) J
u (t)
(t r
dr
.(11)
Аналогично формуле (10), для интервала (0,1) на интервале (1,2) получим модифицированную формулу Герасимова-Капуто (12): 1 <а< 2, п = 2
feu )(t) = -
1 I \d2u(r) dr u (0)
Г(2 -a)[J dr2 (t -r)a-1 + tY1
.(12)
(. (7)
Заметим, что формула (10) переходит в функцию и (г) при а^+0, а формула (12)
du(t)
переходит в производную- при а ^ 1 + 0.
dt
Определение 2. Модифицированной формулой Герасимова-Капуто на интервале 0 <а< 1, п = 1 определим формулой (10), а на интервале 1 <а< 2, п = 2 определим формулой (12) .
Определение 3. Модифицированной формулой Герасимова-Капуто на интервале п-1 <а< п, п = [а]+1 определим формулой
Сравнивая формулы (5) и (7), получим Я! (^ а ^ 0 + 0) = 1, а е (0,1). Также добавим в формулу (3) слагаемое вида
1 ( trdnu{f) dT u("-1)(0) ^ Г(n-a) drn (t-r)a-n+1 + ta
. (13)
( (8)
Формула (13) переходит в формулу и<п при а ^ п -1 + 0.
Рассмотрим обобщение производной дробного порядка от степенной функции:
g2 (t,a^ 1 + 0) = 1, a e (1,2). Преобразуем интеграл в формуле (2) по частям
dkxn dxk
= n(n-1)...(n-k + 1)xn-k = n! xn-k .
(n - k)!
/
0
1
0
Для производной дробного порядка а получим аналогичное выражение:
с1ах"
п!
а _ Г(п + 1) а
Сха (п - а)! Г(п - а + 1)
Разложим гладкую функцию в ряд Мак-лорена, используя предыдущую формулу:
^ и (п)(0)Ьп и(п)(0уп
и(Ь) = Х =Х
п=0
п! п=0 Г(п + 1)'
к,и 4)
п=0 Г(п - а +1) Г(п +1) п=0 Г(п - а +1) ™ (—1)к
Для функции и(Ь) = С08(Ь) = Х-Ь2к получим
к =0 (2к)!
(%ь С08(Т)У
(1()
1
Г (1 - а) к=
да ( 1\к* 2к- а
X , « е (0,1)
П ('■ - а)
1
Г(2 - а) к=1
да ( 1\к* 2к - а
X , а е (1,2)
П ('■ - а)
™ (—1)к
Для функции и (/) = 8т(Ь) = ——-/2к+1 имеем
(о^И^))-
1
Г (1 - а) к=
к=0 (2к +1)
да / 1\к, 2к+1- а
XX , а е (0,1)
П ('■ - а)
1
Г(2 - а) к=0
да ( Х\к + 2к+1- а
X , а е (1,2)
П ('■ - а)
множитель
Г (1 - а )
и слагаемое
I
I (и(Ь)) = [
и (г)Сг _ г и (г)Сг г и (г)Сг (¡-гУ =' (ь-г)а + [ (ь-г)а
В работе [9] получена составная квадратурная интегральная формула с равномерным шагом и с 12 порядком погрешности
о(к12), которую мы используем для первого интеграла в (17) на отрезке [0,Ь]
[и(х)Сх = 5ЬХХС1и(х1)+о(кП),п = 10т,к = Ь——,т е N ,
если
С =
.(15)
16067
299376 16067
1(9688 26575
7(8(4 -16175
99792 5675
-, 1 = 0 VI = п
,(г =0гаэа10)л(0 <г <п) (г = 1гаэа10)V (г = 9гаэа10) (г = 2гаэа10)V (/ = 8гаэа10) . (18)
6237 - (825
55(( 17807
12(7(
(г = 3 газа 10)V (/ = 7 газа 10) , (г = (тза 10) V (г = 6тоа 10) г = 5тоа10
1-1
Вычислим второй интеграл в (17) [
и (Ь - х)Сх
.(16)
Предварительно введем вспомогательный интеграл с параметром a под аргументом функции:
12 (Ь,Ь, а, а ) = |
и(Ь + а - х)Сх
(19)
Опишем численный алгоритм, аппроксимирующий формулу (10), отбрасывая известный 1 и(0)
Представим интеграл (2) в виде суммы, используя во втором интеграле замену переменных х\= Ь - г\'ь, ёх = -ёг.
г и (г)Сг г и (Ь - г'С г и (г)Сг '—Ь и (Ь - ¿)Сг .
=[ ЙГ 1 (17)
Для максимального упрощения алгоритма, первую производную в (17) от известной функции считаем заданной.
Параметр a необходим для аппроксимации первой производной в (10) квадратурной формулой, которая в данной работе не используется, но получена в работе [9].
Весовая функция неотрицательна
р(г) = -1 > 0 на отрезке [0,ь - Ь]. 2
Построим [8] квадратурную формулу Гаусса с тремя узлами.
Утверждение 1. Если функция и(ь) непрерывно дифференцируема на отрезке
[0,1-Ь] с весом р(г) = -1 > 0,г е [0,Ь - Ь], то
г 0
квадратурная формула Гаусса имеет вид
' [и(Ь + а — г)Сг = с^)+ сАхг)+ С}и(х})+ о((ь — Ь)6), (20)
0 2<0
где неизвестные С(,С2,С3, х1, х2, х3 равны
а
1=1
1=2
а
г
1=1
а
2
0
1=2
а
0
С =
(х1 х2 )(х1 Х3)
( 0 - л)1-
(* + а - х)(* + а - х) --
V
(1 -«)
А, л \(*-')2-а (*-'У-аЛ
-(2(! + а) - х - Х Г--— + Л--—
(2 -а) (3 - а)
Л* - ь)1-а
М + а - х)(х - * - а)-----+
31 (1 -а)
С =_:_
(х1 - Х2 )(х2 - Х3 )
ь, ч ^ - ь)2-а (* - ь)
+ (2(* + а) - х - Х Г------—
(2 - а) (3 - а)
3-а Л
-Ь)1-а с _с С3 (1 -а) С2 С1'
„ ь3 (< ь
х = * + а - 2. —соя — | + —, 1 V 3 1 3 I 3
„ Ь, (< + 2я\ Ь,
х = * + а - 2.—соя --| + —,
2 V 3 1 3 I 3
„ Ь3 (< + 4ж) Ь х = * + а - 2, — соЯ --| + ,
3 V 3 \ 3 | 3 '
(
< = агссоя - — | = агссоя
|ь3 - Ьс+а
27 11 3
В формуле (20) Ь1, с, а - неизвестные для ортогонального полинома [8]:
Р (г) = г3 + Ь,г2 + сг + а, р(г) = -1 > 0,7 е [0,* - Ь],
, - 3(3 -а)(* - Ь) 3(2 -а)(3 -а)(* - Ь)2 Ь =-, с =
а =
(6 -а) (5 -а)(6 -а)
- (1 - а)(2 - а)(3 - а)(* - Ь)3 (4 -а)(5 -а)(6 -а)
(21)
Доказательство. Определим ортогональный полином
г3 + Кг2 + сг + а,г|0 = *-Т >0
1 ' 1/-Ь 1Ь
с 3 узлами и р(г) = ^^- > 0,7 е[0, * - Ь] с помощью системы уравнений [8]:
-Ь I-Ь 3,12,, I с с г + Ьг + сг + а | р(г)Р (г)аг = 0 о |-^-аг = 0
Ч а п *Г(г3 + \г2 + сг + а)г ] р(г)?ъ(г)гаг = 0 о ] ^-^-
аг = 0 о
[ р(г) ,( гга = 0 о]' ^^
аг = 0
(* - Ь)4-а + '1 (* - Ь)3-а с(* - + Ь)2-а а и - + Ь)1-а = 0
(4 -а) (3 -а) (2 -а) (1 а)
(* - Ь)5-а + '1 (* - Ь)4-а с(* - + Ь)3-а а с - + Ь)2-а = 0
(5 -а) (4 -а) (3 -а) (2 -а)
(г - Ь)6-а + '1 (* - Ь)5-а с(* - + Ь)4-а а (* - + Ь)3-а = 0
(6 -а) (5 -а) (4 -а) (3 -а)
(24):
Из системы (22) получим уравнения (23),
(* - Ь)3
(4 -а)(2 -а)
(* - ь)3
[(5 - а)(1 - а)
(* - ')3
(5 -а)(3 - а)
(* - Ь)3
, (* - Ь)2 с(* - Ь) + Ь—-—-— +——'Т+-
+ Ь:
1-а)(2-а)
а
-% — /у и — {
- = 0
(* - Ь) с(г - Ь)
, (* - Ь)2 с(* - Ь)
+ь —-—-—+—-—-—+-
- = 0
а
(6 -а)(2 -а) 1 (5 -а)(2 -а) (4 -а)(2 -а) (2 -а)(3 -а)
Вычитая из вторых уравнений систем (23) и (24) первые уравнения в (23) и (24) соответственно, получим систему двух уравнений (25) с константами Ь1, с :
(* - ь)2
(* - 'У
3
(5 -а)(4-а)^
3 ' (6 - а)(5 - а)
+ (* -'
+ (* - ')'
2
| + с
(4-а)(3-а) I 1 (3-а)(2-а)
1
2
| + с
(5-а)(4-а) I I (4-а)(3-а)
1
= 0
(25)
Второе уравнение в (25) умножим на
дробь —1—, затем вычтем из полученного
(2 -а)
выражения с множителем
1
уравнение
3(* -')
(
(5 -а) 2'
1
(4 -а) 1
первое
(6 -а)(2 -а) (4 -а)2
1
1
(4-а) V(5-а)(2-а) (4-а)(3-а)
= 0 о
3(* -')
(5 -а) 2'
16 - 8а + а - (12 - 8а+а") (6 -а)(2 -а)(4 -а)2
2Л
(4 -а)
Г12-7а+а2 -(10-7а + а2 Л (5 -а)(3 -а)(2 -а)(4 -а)
= 0 о
3(* -')
(
4
(5 - а) ^ (6 - а)(2 - а)(4 - а)2
2' (_2_
(4 - а) V (5 - а)(3 - а)(2 - а)(4 - а)
3(* - Ь)(3 -а)
= 0 о
'1 =-:
(6 -а)
(26)
а
1
а
0
+
+
0
2
0
+
+
+
+
0
0
а
Выражая из первой формулы системы (25), найдем с с учетом результата (26):
с = -(3 -а)(2 -а)1 ^ - Ь)
-3-| + (t - Ь)Ь,|-2-
(5-а)(4-а) I I (4-а)(3-а)
3(3 -а)(2 -а)(1 - Ь)2
1
2
(4-а) ^ (5-а) (6-а)
3(3 -а)(2 -а)(t - Ь)2 / ч 3(2 -а)(3 -а)(t - Ь)2
(4-а)(5 -а)(6 -а)
■(- 4-a)—-
(5-а)(6-а)
. (27)
3р - Ь)2(3 -а) (6 -а)2 (5 -а)
(15 - 8а + а - (12 - 8а + а2)) =
9^ - Ь)2(3 -а) (6 -а)2 (5 -а)
> 0 :
так как 0 < а < 1.
Поэтому воспользуемся формулами А.Д. Фаддеева [7, стр. 65] для трех действительных различных корней уравнения (29). После-довательно вычисляем:
г, р, у, у2, у3, , , , х1, х2, х3 г =
Получим й из первой строки системы (22) с учетом (26), (27) для найденных Ь1, с :
d = -(1 -а) = -(1 -а^ - Ь): - (1 -а)(t -Ь)3
Г(( -Ь)3 | ъ ^ -Ь)2 | С(( -Ь) (4-а) 1 (3-а) (2-а)
1 3 3(3-а)
Л
- + -
(4-а) (6-а) (5-а)(6-а) 6-а-(12 - 3а) 3(3-а)
(4 -а)(6 -а) (5 -а)(6 -а)
2
3
_ (1 -а)(3 -а)(t - Ь)3 = (6-а) ^ (4-а) (5-а)
(1 - а)(3 -а)^ - Ь)3 (10 - 2а - (12 - 3а)
(6-а) ^ (4 -а)(5-а)
(1 -а)(2 -а)(3 -а)(t - Ь)3 (4 -а)(5 -а)(6 -а) '
(28)
Получен ортогональный полином с
учетом формул (26), (27), (28):
Р ф = 23 + Ь22 + + d = 23 -
3^ - Ь)(3 -а) (6 -а)
2 -
3(2 -а)(3 -а)(t - Ь)2 (1 -а)(2 -а)(3 -а)(t - Ь)3 псл +-;-;-2--;-;-;- . (29)
(5 -а)(6 -а)
(4 -а)(5 -а)(6 -а)
р = агссо^ -~\ = агсс°9
2 V 3 Ь{
Ь,' —^ + d 27 I 1 3
(30)
~ Ь3 (р| „ Ь3 ( р + 2п '
У3=^Ч^). 2.=^Я?)-*
Ь (р-Ь Ь (р+ Ь
2 = 2.1— соя ■
2 V 3 \ 3
А также
■
—Ч 2 = 2Л — соя -
3 1 3
Поскольку корни ортогонального полинома действительны [8], попарно различны, то уравнение (29) имеет три положительных корня на отрезке [0^-Ь]. Кубическое уравнение (29) согласно работе [7] следует привести к каноническому уравнению заменой переменных:
Ь1 ______л „ Ь12 ( 2 \ 3 Ь1С
X = t - а - ^ = t - а - 2^-Ь1-соя(| + у,
„ Ь,3 (р + 2жЛ Ь, , ч X = t - а - ^ = t - а - 2у-^соа —)-у. (31)
„ ,Ь, (р + 4ж\ Ь, X = t - а - ^ = t - а - соя ^—) + у
Весовые коэффициенты С, С, С в формуле Гаусса с найденными корнями (31) найдем с помощью классической задачи [8]:
ч , u(t + а - z)dz t- dz (- Ь)1-" „ „ „
+ а - 2) = 1: Г -(-аГ2±~ = I — = „ \ = С + С2 + С3
0 2а 0 2а (1 -а)
г (t + а - „ г dz г
+ а - 2) = t + а - 2 : Г --^— = (t + а) Г — - Г а
^(t - Ь) а ^ - ЬТа -г +Г
— (t I а) — I с2^2 +
(1 -а) (2 -а)
у = 2-у,у -ру-Ч = 0,р = Сч = ^)Ь1 -у-d.
По критерию Д.К. Фаддеева [7] кубическое уравнение (29) имеет три раз-личных
^ 2
вещественных корня если р1 = -р = ^^— с > 0 .
Это условие выполнено и для нашего ортогонального полинома (29).
_Ь12 _3(/-Ь)2(3-а) ( (3-а) (2-а)
Р1 =--с =-1---
3 (6-а) I (6-а) (5-а)
\2 (- а - 2)2 <¡2 2 С ¿2 ' zdz 22<!2
—(- а -1) = ( - а - 2) : | ^-^-= (/ - а)2 I —-2(/ - а) I —-
2 с 2 с 2 с 2
- а - х) = (( - а - 2) : \ -
0 2 0 2 0 2 0 \2-а Л »\3-а
- а)2 ((г-)1а-2(t - а) — С1х12 - С2 х22 - СЛ2 О
(1 - а) (2 - а) (3 - а)
' (I-Ь }-а_
(1 -а) = С'- С2 - С3 л((-Ь)1-а к-Ь)2-а
t - (Та)" ^0)= СХ - С2 Х2- Сз Хз
, ,2 ^-ЬМ--)2гa ^-ЬТ" г г г г Г г П - а)^—'-— 2и- ау—'— ^—'— = Сх - С,х,2 - С,х, (1 -а) (2-а) (3-а) 1 1 2 2 3 3
'.(3 2)
3
3
2
3
0
Первую строку в системе уравнений (32) умножим на х3 и вычтем из второй строки:
(' + а - х3)- Ь) = C1 (х1 - х3)+ С2(х2 - х3) . (33)
(1-а) (2-а) Затем вторую строку системы (32) умножим на число х3 и вычтем из третьей строки:
, 42 (' - Ь)1-" М— ь)2-- ('- ь)3-" ч (' - ь)1-- х3 ('- ь)2-0
(' + а)-----2(1 + ан--— + ^----(' + а)х, ---— + —-—
(1 —а) (2-а) (3-а) 3 (1-а) (2-а)
= Сх (х - х)+С2х2 (х - х )о (I+а)('+а - х)
(' - Ь)1-
-(2(*+а) - х
(2-а) (3-а)
(1 -а)
=Сх (х - х3)+ Сгхг (х - х )((34)
Выражение (33) умножим на х2 и вычтем из полученного (34).
, V л' - ')1—: ^ ч ч(/ - Ь)2-: (' - ь)3-:
(! + а)(' + а - х )Л----(2(г + а) - х )----+ Л----
3 (1 -:) у ' (2-:) (3-:)
_ ч (' - ь)1-: (' - ь)2-: _
(' I а х3 )х^ + х^ —
(1 -а) (2 -а)
= С1х1 (х1 — х3 )— С1х2 (х1 — х3 )= С1 (х1 — х2 )(х1 — х3 ) =
Л' - ¿У ^ ^ - Ь)2-" (' - Т
= (I + а - х )(' + а - х )-----(2(1 + а) - х - х Г-^—+ Л--— ■'
(1 -а) (2-а) (3-а)
Из последнего уравнения выразим С :
1 (, ч(' - Ь)1—" (' + а - х)(' + а - х) ---
С, =
(х1 х2 )(х1 х3)
(,, ч ^—ь)2-а ('—ь)3-аА
— (2(1 + а) — х — X ) - — +-
(1 -а)
. (35)
(2 — а) (3 — а) С учетом (35) из формулы (33) найдем С2 :
С1 (х1 - х3 ^ 1
£ _ С1 (х1 хд/ +_
(х2 - х3) (х2 - х3)
( (' - ь)1-а (' - ь)2-а ^
(' + а - х р—---А—-
(I + а - х)
(х2 х3)
(I+а - х )(х -' - а) (' - Ь)-" (2('+а) - х, - х)(' - Ь
(1 -а) (2-а)
' - 'Т-
(2 -а)
(' - ь)-а (' - ь)2-аЛ
(1 -а) (2 -а)
--.-1-, I (' + а - х2)(( + а - х3)' - ' -(2(( + а) - х3 - х2 ^ ' + ' ' | +
(х1 -х21х2 -х3(1-а) (2-а) (3-а) |
(' - Ь
(х - х )(х - х) (1—-) (х - X )(х - х) (2-") (х - х )(х - х )(3-")
(
1
(х1 х2 )(х2 х3)
\1-:
— и I —
('+а - х3)(х1-'- +(2('+- х- х)- , „ ,
(1-:) (2-:) (3-:)
Из системы (32) выразимС :
_мт (1 -:)
С3 =-
—С — С
С2 С1
(36)
(37)
Результаты формул (37), (36), (35), (31), (30), (29), (28), (27) совпадают с формулами (20), (21) и Утверждение 1 доказано. Отметим также, что в упрощенном алгоритме квадратурная формула Гаусса (20), (21) для вычисления второго интеграла в формуле (17)
(или в формуле (19) с параметром a=0) используется одни раз. Как и используется один раз формула (18) с равномерным шагом для вычисления первого интеграла в (17).
Замечание 1. Для дробной производной на интервале (0,1) в упрощенном алгоритме в Утверждении 1 нужно выбрать функцию и вес: и( г - г) ^ и (г- г ),а = 0,р(г) = 1/г:,а е (0,1).
А на интервале (1,2) нужно выбрать функцию и (х) с другим весом:
и(г - г) ^ и (г- г), а = 0, р(г) = 1/ г:-1, а е (1,2) .
И в этом случае две формулы (20), (21) Утверждения 1 останутся верными, если в них формально заменить а-->а-1 и заменить весовую функцию.
Пользуясь формулами (10), (12) и соответствующим численным алгоритмом (18), (20), (21), составим таблицу для численного значения дробной производной и аналитического значения производной с помощью ряда (16) для функции / (I) = ят('), I = 1,: е (0,2). Таблица 1. Производная Герасимова-Капуто для функции и(г) = ятО,' = ж/2,ае (0,2)
а (о:+,,яш(о) (д:+,ят(г)) \ и+,' \ / >пит л(о0:+,,)
(| = я-/2) (| = я72) 'оТ п
10"14 0.9999999999 0.99999999999 1.7763568
99995 9994 3Е-015
10"8 0.9999999952 0.99999999527 -4.1078Е-
79993 9997 015
0.1 0.9470275154 0.94702751540 5.9952043
07809 7803 32Е-015
0.4 0.7184255171 0.71842551711 1.4432899
19136 9135 32Е-015
0.5 0.6197920300 0.61979203007 -1.4432Е-
73509 3510 015
0.6 0.5109932063 0.51099320638 2.2204460
84873 4873 49Е-016
0.9 0.1357686082 0.13576860826 2.4980018
62797 2797 05Е-016
1-10"8 1.3707621751 1.37076218024 -5.0526Е-
9046Е-008 3091Е-008 017
1-10"14 1.3871023850 1.37579260713 1.1309777
6829Е-014 4274Е-014 93Е-016
1+10- -1.371991 Е- -1.3100631 Е- -6.1927Е-
14 014 014 016
1+10-8 -1.370762Е- -1.3707621Е- -3.7398Е-
008 008 016
1.1 -0.13758466 -0.13758466 3.9690473
5279301 5279305 13Е-015
1.4 -0.53501308 -0.535013081 -2.2204Е-
1857711 857709 015
1.5 -0.65277766 -0.65277766 -8.8817Е-
5939620 5939619 016
1.6 -0.75817951 -0.758179518 -6.6613Е-
8757917 757916 016
1.9 -0.97151962 -0.971519622 -2.3314Е-
2786878 786876 015
2-10-8 -0.99999999 -0.999999998 - 2.220Е-
8353808 353809 016
2-10-14 -0.99999999 -0.999999999 1.1102230
9999998 999998 2Е-016
а
1
2:
3-а
Аналогичную таблицу для численного и аналитического значений производной с помощью ряда (15) составим для функции f (t) = cos(t), t = 1,ae (0,2). Первая часть и вторая часть каждой таблицы на интервалах (0,1) и (1,2) находится разными алгоритмами и программами с особенностями, описанными в Замечании 1.
Таблица 2. Производная Герасимова-Капуто для функции u(t) = cos(t), t = г/2,ае (0,2)
a fe cos(t)) Da+t coso) A(D0a+, cos(t
(t = —/2) (t = —/2) (t = —/2)
10-14 -1.36710160 -1.310063169 -5.70384E-
44627E-014 057677E-014 016
10"8 -1.37076214 -1.370762118 -3.18475E-
99648E-008 117236E-008 016
0.1 -0.13758466 -0.137584665 4.218847493
5279301 279305 57E-015
0.4 -0.53501308 -0.535013081 -1.77635E-
1857711 857709 015
0.5 -0.65277766 -0.652777665 -8.88178E-
5939620 939619 016
0.6 -0.75817951 -0.758179518 -1.22124E-
8757917 757916 015
0.9 -0.97151962 -0.971519622 -2.22044E-
2786878 786876 015
1- -0.99999999 -0.999999998 2.220446049
10-8 8353808 353809 25E-016
1-10- -0.99999999 -0.99999999 1.110223024
14 9999998 9999998 625E-016
1 + -0.99999999 -0.9999999 1.110223024
10-14 9999998 99999998 625E-016
1 -0.99999999 -0.999999995 3.885780586
+10-8 95279994 279997 188E-015
1.1 -0.94702751 -0.947027515 -5.99520E-
5407809 407803 015
1.4 -0.71842551 -0.718425517 -1.33226E-
7119136 119135 015
1.5 -0.61979203 -0.619792030 1.443289932
0073509 073510 012E-015
1.6 -0.51099320 -0.510993206 -5.55111E-
6384873 384873 016
1.9 -0.13576860 -0.135768608 -3.33066E-
8262798 262797 016
2- -1.37076217 -1.370762165 -1.01881E-
10-8 51904E-008 002316E-008 016
2- -1.38710238 -1.375792607 -1.13097E-
10-14 50682E-014 134274E-014 016
Таблицы 1 и 2 получены при следующих параметрах п=10000^=Ь (п-число интервалов для отрезка в первом интеграле в (17)), !=45, 45h=t-b -длина правого отрезка в квадратуре Гаусса.
Первый столбец - порядок производной а, второй столбец - аналитическое значение производной Герасимова по формулам (15), (16), третий - численное значение производной по алгоритму (18), (20), (21). Четвертый столбец - разность между аналитическим и численным значениями дробной производной.
Из табл. 1 для функции u(t) = sin(t),t = —/2,ae[0,2],sinf—| ^ 1 -0
V 2 ) a^+0
cosí — I ^ + 0,-sin — I ^ -1 + 0,
V 2 )a^1-0 V 2 )a^2-0
что соответствует производным u(t) = sin(t). Из табл. 2 для функции
u(t) = cos (t), t = —/2,ae[0,2],cosf—J = -0 ^a =+0 -sinf — I = -1 ^a = 1,-cos — I = -0^a = 2-0 ,
V 2 J \ 2 J
что соответствует производным u(t) = cos(t). Основные полученные результаты:
1) Впервые получены модифицированные формулы (10), (12), (13) Герасимова-Капуто.
2) С учетом модифицированных формул (10), (12) получен численный алгоритм (18), (20), (21) и доказана корректность алгоритма - Утверждение 1.
3) Результаты алгоритма (18)-(21) сравнены с аналитическими формулами (15), (16), невязка вычислений не превышает 10-15 .
4) Формулы (20), (21) Утверждения 1 остаются верными при формальной замене в них a на интервале (0,1) на а-1 на интервале (1,2) .
Список источников
1. Корчагина А.Н. Использование производных дробного порядка для решения задач механики сплошных сред // Известия Алтайского государственного университета. 2014. № 1-1(81). С. 65-67. DOI 10.14258/izvasu(2014)1.1-14. EDN SECUCD.
2. Нахушев А.М. Дробное исчисление и его применение. М.: Физматлит, 2003. 272 с.
3. Бештокова З.В. Устойчивость и сходимость
монотонных разностных схем, аппроксимирующих краевые задачи для интегро-дифференциального уравнения с дробной по времени производной и оператором Бесселя / З.В. Бештокова, М.Х. Бештоков // Дифференциальные уравнения и процессы управления. 2021. № 3. С. 26-50. EDN GMBWPR.
4. Бештоков М.Х. Краевые задачи для нагруженного модифицированного уравне-ния влагопереноса дробного порядка с оператором Бесселя и разностные методы их решения // Вестник Удмуртского универси-
)
тета. Математика. Механика. Компьютерные науки. 2020. Т. 30, J№ 2. С. 158-175. DOI 10.35634/vm200202. DNHMCSFN.
5. Алиева С.Т., Мансимов К.Б. Условие оптимальности типа принципа максимума Понтрягина в задаче управления линейными разностными уравнениями дробного порядка // Вестник Пермского университета. Математика. Механика. Информатика. 2023. Вып. 4(63). С. 5-11. DOI 10.17072/ 1993-0550-2023-4-5-11. EDN ACKUPX.
6. Мансимов К.Б., Ахмедова Ж.Б. Аналог принципа максимума Понтрягина в задаче оптимального управления системой дифференциальных уравнений с дробной производной Капуто и многоточечным критерием качества // Вестник Пермского университета. Математика. Механика. Информатика. 2022. Вып. 3(58). С. 5-10. DOI 10.17072/1993-0550-2022-3-5-10. EDN THSSNA.
7. Фаддеев Д.К. Лекции по алгебре: учеб. пособие для вузов. М.: Наука. Физматлит, 1984. 416 с.
8. Бахвалов Н.С., Лапин А.В., Чижонков Е.В. Численные методы в задачах и упражнениях. М.: БИНОМ. Лаборатория знаний, 2010. 240 с.
9. Пастухов Д.Ф., Пастухов Ю.Ф., Волосова Н.К., Волосов К.А., Волосова А.К. Вычисление производных дробного порядка с высокой степенью точности. Новополоцк: ПГУ, 2020. 21 с. URL: hhtps://elib.psu.by/ handle/123456789/25335.
10. Гербер А.Д. Описание алгоритма приближенного вычисления несобственного интеграла, определяющего значения дробной производной // Математика и естественные науки. Теория и практика: межвуз. сб. науч. тр. Т. Вып. 16. Ярославль: Ярослав. гос. техн. ун-т. 2021. С. 22-31. EDN CYCCAJ.
References
1. Korchagina, A.N. (2014), "The use of fractional derivatives to solve problems in continuum mechanics", Proceedings of the Altai State University, no. 1-1(81), pp. 65-67. DOI 10.14258/izvasu(2014)1.1-14. EDN SECUCD.
2. Nakhushev, A.M. (2003), "Fractional calculus and its application",M., FIZMATLIT, 272 p.
3. Beshtokova, Z., Beshtokov, M.Kh. (2021), "Vstability and convergence of monotone difference schemes that approximate boundary value problems for an integro-differential equation with a time-fractional derivative and the Bessel operator / Z.V. Beshtokova", Differential equations and control processes, no. 3, pp. 26-50. EDN GMBWPR.
4. Beshtokov, M.Kh. (2020), "Boundary value problems for the loaded modified fractional order moisture transfer equation with the Bessel operator and difference methods for their solution", Bulletin of the Udmurt University. Mathematics. Mechanics. Computer Science, vol. 30, no. 2, pp. 158-175. DOI 10.35634/vm200202. EDN HMCSFN.
5. Alieva, S., Mansimov, K.B. (2023), "Toptima-lity condition of the type of Pontryagin's maximum principle in the problem of controlling linear difference equations of fractional order". Bulletin of Perm University. Mathematics. Mechanics. Computer Science, no. 4(63), pp. 5-11. DOI 10.17072/1993-05502023-4-5-11. EDN ACKUPX.
6. Mansimov, K.B., Bakhmedova, Z.H. (2022), "An analogue of the Pontryagin maximum principle in the problem of optimal control of a system of differential equations with the Caputo fractional derivative and a multipoint quality criterion", Bulletin of Perm University. Mathematics. Mechanics. Computer Science, no. 3(58), pp. 5-10. DOI 10.17072/1993-05502022-3-5-10. EDN THSSNA.
7. Faddeev, D.K. (1984), "Lectures on algebra: Textbook for universities", M., Science, Fizmatlit, 416 p.
8. Bakhvalov, N.S., Lapin, A.V., Chizhonkov, E.V. (2010), "Numerical methods in problems and exercises", M., BINOM, knowledge laboratory, 240 p.
9. Pastukhov, D.F., Pastukhov, Y.F., Volosova, N.K., Volosov, K.A., Volosova, A.K. (2020), "Calculation of production of the debrief order with hign accuracy", Novopolotsk, PGU, 21 p. URL: hhtps://elib.psu.by/handle/123456789 /25335.
10. Gerber, A.D. (2021), "Description of the algorithm for approximate calculation of an improper integral that determines the values of the fractional derivative", Mathematics and natural sciences. theory and practice: Interuniversity collection of scientific papers, Yaroslavl, Yaroslavl State Technical University, vol. 16, pp. 22-31. EDN CYCCAJ.
Информация об авторах:
Наталья Константиновна Волосова - аспирант МГТУ им. Н. Э. Баумана (105005, Россия, г. Москва, 2-я Бауманская ул., д. 5, стр. 1), [email protected];
Константин Александрович Волосов - доктор физико-математических наук, профессор кафедры прикладной математики Российского университета транспорта (127994, ГСП-4, Россия, г. Москва, ул. Образцова, д. 9, стр. 9), [email protected], AuthorlD: 128228;
Александра Константиновна Волосова - кандидат физико-математических наук, начальник аналитического отдела ООО "Трамплин" Российского университета транспорта (127994, ГСП-4, Россия, г. Москва, ул. Образцова, д. 9, стр. 9), [email protected], AuthorlD: 607500;
Михаил Иванович Карлов - кандидат физико-математических наук, доцент кафедры высшей математики Московского физико-технического университета (МФТИ) (141701, Россия, Московская область, г. Долгопрудный, Институтский пер., 9.), [email protected], AuthorlD: 14680;
Дмитрий Феликсович Пастухов - кандидат физико-математических наук, доцент кафедры технологий программирования Полоцкого государственного университета (211440, Республика Беларусь, Витебская обл., г. Новополоцк, ул. Блохина, 29), [email protected], AuthorlD: 405101;
Юрий Феликсович Пастухов - кандидат физико-математических наук, доцент кафедры технологий программирования Полоцкого государственного университета (211440, Республика Беларусь, Витебская обл., г. Новополоцк, ул. Блохина, 29), [email protected], AuthorlD: 405109.
Information about the authors:
Natalya K. Volosova - Post-graduate Student of Bauman Moscow State Technical University (2nd Baumanskaya St., 5-1, Moscow, Russia, 105005), [email protected];
Konstantin A. Volosov - Doctor of Physical and Mathematical Sciences, Professor of the Department of Applied Mathematics of the Russian University of Transport (Obraztsova St., 9-9, Moscow, GSP-4, Russia, 127994), [email protected], AuthorlD: 128228;
Aleksandra K. Volosova - Candidate of Physical and Mathematical Sciences, Chief Analytical Department "Tramplin" LLC, Russian University of Transport (Obraztsova St., 9-9, Moscow, GSP-4, Russia, 127994), [email protected], AuthorlD: 607500;
Mikhail I. Karlov - Candidate of Physical and Mathematical Sciences, Associate Professor of the Department of Higher Mathematics, Moscow University of Physics and Technology (Institutskiy per., 9, Dolgoprudny, Moscow region, Russia, 141701), [email protected], AuthorlD: 14680;
Dmitriy F. Pastukhov - Candidate of Physical and Mathematical Sciences, Associate Professor of Polotsk State University (Blokhin St., 29, Novopolotsk, Vitebsk Region, Republic of Belarus, 211440), [email protected], AuthorlD: 405101;
Yuriy F. Pastukhov - Candidate of Physical and Mathematical Sciences, Associate Professor of Polotsk State University (Blokhin St., 29, Novopolotsk, Vitebsk Region, Republic of Belarus, 211440), [email protected], AuthorlD: 405109.
