Вестн. Сам. гос. техн. ун-та. Сер. Физ.-мат. науки. 2014. № 1 (34). С. 125—133
Математическое моделирование
УДК 519.642.2
АНАЛИЗ РАЗНОСТНОЙ СХЕМЫ АНАЛОГА ВОЛНОВОГО УРАВНЕНИЯ С ОПЕРАТОРОМ ДРОБНОГО ДИФФЕРЕНЦИРОВАНИЯ
В. Д. Бейбалаев1, А. З. Якубов2
1 Дагестанский государственный университет,
Россия, 367025, Махачкала, ул. М. Гаджиева, 43a.
2 Дагестанский государственный институт народного хозяйства,
Россия, 367008, Махачкала, ул. Д. Атаева, 5.
Проведён анализ разностной схемы краевой задачи для аналога волнового уравнения. Исследованы явные и неявные разностные схемы для численного решения первой краевой задачи аналога волнового уравнения с оператором дробного дифференцирования Капуто. Методом гармоник Фурье доказаны критерии устойчивости этих разностных схем. Получены оценки собственных значений оператора перехода с одного временного слоя на другой. На примере проведён вычислительный эксперимент по анализу предложенной разностной схемы. Построены графики численного решения краевой задачи для волнового уравнения с оператором дробного дифференцирования при различных значениях параметров дробного дифференцирования а и в. Установлено изменение периода колебаний при переходе к дробной производной. На примере показано, что параметры а и в становятся управляющими.
Ключевые слова: волновое уравнение, частная дробная производная Капуто, разностная схема.
Введение. Развитие науки и техники привело к широкому использованию на практике самых различных волновых процессов. В результате этого живые системы подвергаются воздействию волн не только естественного, но и искусственного происхождения. Можно сказать, что живые организмы везде сталкиваются с различного типа волнами. Физические поля самых разных типов в процессе эволюции живых организмов оказывают на них существенное влияние. Сравнительно давно появилось отчетливое понимание того, что живые организмы не только активно используют волновые процессы в своей деятельности, но и способны чутко реагировать на внешние воздействия, имеющие волновую природу.
ISSN: 2310-7081 (online), 1991-8615 (print); doi: http://dx.doi.org/10.14498/vsgtu1209 © 2014 Самарский государственный технический университет.
Образец цитирования: В. Д. Бейбалаев, А. З. Якубов, “Анализ разностной схемы аналога волнового уравнения с оператором дробного дифференцирования” // Вестн. Сам. гос. техн. ун-та. Сер. Физ.-мат. науки, 2014. № 1(34). С. 125-133. doi: 10.14498/vsgtu1209. Сведения об авторах: Ветлугин Джабраилович Бейбалаев (к.ф.-м.н., доц.), доцент, каф. прикладной математики. Амучи Загирович Якубов (к.ф.-м.н., доц.), доцент, каф. информатики и вычислительной техники.
E-mail addresses: [email protected] (V.D. Beybalaev, Corresponding author), [email protected] (A.Z. Yakubov)
125
В. Д. Бейбалаев, А. З. Якубов
Ответная реакция живой системы на внешние волновые воздействия может происходить на различных структурных уровнях живого организма — от молекулярного до субклеточного, клеточного уровня. Волновые процессы активно используются и в технической деятельности людей.
Несмотря на значительные усилия исследователей, до сих пор задача создания адекватных количественных моделей волновых процессов остаётся актуальной. Особенно это актуально, когда речь идет о системах с фрактальной структурой. При описании свойств систем с фрактальной структурой нельзя использовать представления евклидовой геометрии, поэтому необходимо привлечь представления геометрии дробной размерности. Особенность систем с фрактальной структурой в том, что для них существенны такие эффекты, как память, сложная природа пространственных корреляций и эффекты самоорганизации. Создание адекватных математических моделей для систем, где проявляются свойства самоорганизации, детерминированного хаоса, также требует привлечения нетрадиционных подходов, основанных на применении математического аппарата дифференциальных уравнений дробного порядка [1-7].
В качестве математической модели волнового процесса в нелокальных средах исследуем в области
D = {(x,t) :0 ^ x ^ l, 0 ^ t ^ T}
первую краевую задачу для волнового уравнения с производными дробного порядка.
Задача. Найти решение уравнения
cD‘0tu(x,t) = C (x, t)cD0x u(x,t) + f (x,t), (1)
удовлетворяющее начальным условиям
u(x, 0) = pi(x), ut(x, 0) = ^2(x) (2)
и граничным условиям
u(0,t) = Pi (t), u(l,t) = P2 (t). (3)
Здесь 1 < a ^ 2, 1 < в ^ 2, f (x,t) — достаточно гладкая функция, C(x,t) ^ 0,
CDotu(x, t)
1
Г(2 - a)
t
0
u"(x, s) (t — s)a-i
ds
— частная дробная производная Капуто [10].
Численный метод решения краевой задачи (1)-(3). Для численного решения задачи (1)-(3) в области D введём сетку с шагом h по x и т по t:
етhT = {(xm,tn) : xm = mh, Ьп = пт; h = l/M, т = T/N;
m = 0,1,... ,M, n = 0,1,... ,n\.
Для дробных производных Капуто в случае 1 < а ^ 2, 1 < в ^ 2 на частичных отрезках [tn-i,tn], [xm-i,xm] имеют место разностные аппроксимации [11]
126
Анализ разностной схемы аналога волнового уравнения ...
1
Г(2 - а)
t
n
tn- 1
u"(x, s) (tn - s)a~l
ds =
u(x, tn+l) - 2u(x, tn) + u(x, tn-1) Г(3 — а)та
+ O (т4-а), (4)
1 fXm u"(s,t)
Г(2 - e) Jxm-1 (xm - s)e-1
ds =
u(xm+1, t) 2u(Xm,t) + u(xm-1, t) t c\f и4-в\ /r. \
Г(3 - e)he +O(h )- (5)
Обозначим через um ~ u(xm,tn) приближённое решение в точке (xm,tn), n = 0,1,... ,N, m = 0,1, . . . , M, аналогично обозначим Cm = C (xm,tn), fm = = f (xm, tn). Воспользовавшись (4) и (5), для уравнения (1) получим разностную схему с весами
u
n+1_2un m - m
+u
n- 1
m
Г(3 - а)та
аЛв um+1 + (1
- а)лвum + f
n,
m,
(6)
где
un+1 — 2un+1 + Iin+1 un — 2un + un
двun+1 = cn um+1 2um + um-1 двun = cn um+1 2um + um-1
Л um = Cm Г(3 - в)he , Л um = Cm Г(3 - e)he ■
В случае a = 0 получаем явную разностную схему на 5-точечном шаблоне:
u
n+1_2un + un-1
m - m m
Г(3 - а)тс
= Cmn
un 2 un + un
um+1 2um + um-1
1 Г(3 - e)he
+ f
(7)
Теорема 1. Явная разностная схема (7) устойчива, если СГ(3 - а)т°
Г(3 - в)he
< 1,
где С = max |С (x,t)l.
Доказательство. Решение un+1 можно представить в виде un+1 = Aun + Г(3 - а)тafn.
Тогда S = A — оператор перехода с одного временного слоя на другой. Условием устойчивости по начальным данным разностной схемы (6) является ||S|| ^ 1. Действительно,
Иуп+Ч = l|Ayn|l < l|A||-||yn|l
или
l|yn+1| = llyn| < ... < lly0|.
Следовательно, начальные возмущения затухают. Условие ||S| ^ 1 есть условие того, что спектр оператора S лежит внутри круга единичного радиуса на комплексной плоскости. Это означает, что max |Л| < 1, где Л — собственное число оператора перехода S.
127
В. Д. Бейбалаев, А. З. Якубов
Для нахождения собственных значений оператора перехода решение um представим в виде возмущения:
um = \ne
ПЛ1ЛШ
где i — мнимая единица, Л — собственные числа оператора перехода [12]. В результате подстановки (8) в (7) получим квадратное уравнение
Л2 — 2(1-
2Cmг(3 - а)т“ о ш
Г(3 - в)ЬР ~~ 2
sin2 Л + 1 = 0.
(9)
Произведение корней уравнения (9) равно единице, а его дискриминант
d = "• sin^t „у; sin^ - Л. (10)
16Cm Г(3 - а)та . 2 Ш ( cm Г(3 - а)та . о Ш_ Л Г(3 - e)he Sin 2V Г(3 - e)he Sin 2 1.
Если дискриминант (10) отрицательный, то корни Л1, Л2 — комплексно-сопряжённые и равны единице по модулю. Следовательно, условие |Л| < 1 выполняется, если они выбраны так, что
СГ(3 - а)та Г(3 - в)he
< 1.
□
Итак, решение um разностной задачи
un+1 _ 2un + un 1 un _ _ 2un + un .
um 2um + um cn um+1 2um + um-1
cm
+ fn
m
Г(3 - а)та ~m Г(3 - e)he
um = ^i(xm), um = um + ^(^m), un = ^1(tn), uM = +2 (tn),
можно найти явным методом.
В случае a = 1 получим полностью неявную схему
u
n+1 _ 2un + un-1 m - 2um + um
m + um __ cn um+1
Г(3 - а)т
n+1 _ 2un+1 + un+1 m+1 " ^m ""г ^m— 1
Г(3 - e)he
+ fn
m
(11)
Теорема 2. Неявная разностная схема (11) безусловно устойчива. Доказательство. Для нахождения собственных значений оператора перехода решение u'm представим в виде возмущения (8). В результате подстановки (8) в (11) получим квадратное уравнение
[1 + р] ■ Л2 - 2Л + 1 = 0,
где
Р
4CmГ(3 - а)та Г(3 - e)he
sin2
Ш 2.
Дискриминант этого квадратного уравнения отрицательный. Таким образом, корни комплексно-сопряжённые и имеет место неравенство
|Л|
Уу
1 + р
< 1.
128
Анализ разностной схемы аналога волнового уравнения ...
Следовательно, неявная разностная схема (11) безусловно устойчива. □ Для n = 1, 2,..., N имеем
(С D0>(x,t))n =
1
ftn „.//
Г(2 - a) J0 (tn - s)
u"(x, s)
a— 1
ds =
nt
t
I ■ u"(x,s)
Г(2 - a) i~1-lti-1 (tn - s)a-1
ds =
Г(2 - a)
ds
1 (u(x,ti+1) - 2u(x,ti) + u(x,ti-i) 2л\ A
— ^----------------------T2--------------- + °<T )) У (tn - s)«-1
i=1
n
1 y fu(x, ti+i) - 2u(x, ti) + u(x, ti-i) +O(T2)ЛГ (tn - s)2 alti
r(2 - a) 1
2 - a Jti_1
Г(3 - a)
Eu(x ti+1) 2u(x, ti) + u(x ti-1) ~ . 2\ A
----------------2--------------- + 0(т ^ X
i=1
X ((tn - ti-1)2 a - (tn - ti)2 a) =
Г(3 - a)r2 ^
y [u(x, ti+1) - 2u(x, ti) + u(x, ti-1)] (tn-i+1- tn-i)+
i=1
+ o(t2)£>n-+1 - tn-) =
Г(3 - a)r2
^ [u(x, tfc+1) - 2u(x, tk) + u(x, tfc-1) x
fc=1
(tn-k+1 - tn-k) + O(T3-a), (12)
1
1
так как tn-i+1- tn-i <т
Аналогично для m =
2- .
1, 2,..., M имеем
(C Dteu(x,t))m
1
Г(2 - в)Уо 1
Г(2 - в)
X
m
u"(s, t) (xm - s)e-1
ds =
u"(s, t) (xm - s)e-1
ds =
1
r(3 - e)h2
[u(xk+1, t) - 2u(xk, t) + u(xk-1, t)] x
k=1
Xx
2-в
m-k+1
-x
2-в
m- k
+ O(h3-e), (13)
Воспользовавшись разностными аппроксимациями (12) и (13), получим составную разностную схему с весами:
1
Г(3 - a)r2
n
ХХ^1
k=1
2u
k
m
+u
k-1
m
t
2-
n- k+1
t2- =
tn k =
129
В. Д. Бейбалаев, А. З. Якубов
= a • Лв+ (1 - а)Лвunm + f
П
m,
где
Лв um+1
C ГУ
1
Г(3 - e)h2
m
£(<+!
k=1
— 2 u
n+1
k
+ u
n+1
k-1
)(x
2-в
m- k+1
— x
2-в
m- k
5
Лв un =
m
Г(3 - да2 k=1
(uk+1 2uk + uk-1 {xm-k+1 Xm—
1
Пример. Найти решение задачи
CD0atu(x,t) = CD0exu(x,t), 0 < x < 4, 0 < t < 10,
u(x, 0) = sin(nx/4), ut(x, 0) = 0, 0 ^ x ^ 4, u(0, t) = 0, u(4,t)=0, 0 ^ t ^ 10.
(14)
На рис. 1-3 приведены графики решений задачи (14) при различных значениях параметров а, в.
Как видно из рис. 1, при переходе к дробной производной по времени период колебаний становится меньше, а при переходе к дробной производной по координате, наоборот, период колебаний растёт. При этом длина волны остаётся неизменной (рис. 2 и 3). Следовательно, при переходе к дробным производным в первом случае скорость распространения волны становится больше, чем в стационарном случае, а во втором случае — меньше.
Выводы. Значительные успехи при рассмотрении волновых процессов в системах с фрактальной структурой связаны с использованием формализма интегралов и производных дробного порядка. Повышенный интерес к дифференциальным уравнениям дробного порядка обусловлен их физической интерпретацией [7,13]. Было показано, что переход к производной дробного порядка по времени позволяет учитывать эффекты памяти системы [7]. Это и вызвало их широкое применение в практически во всех областях естествознания.
Область применимости решений дифференциальных уравнений дробного порядка значительно шире, чем дифференциальных уравнений с целочисленным дифференцированием, поскольку последние оказываются их частным случаем. Уравнения в производных дробного порядка позволяют учесть обратимые и необратимые процессы [7,13]. Это позволяет не только получить принципиально новые, но более глубоко осмыслить известные результаты, позволяя при этом создать адекватные количественные модели исследуемых явлений
Исследование свойств систем с фрактальной структурой должно быть проведено с учётом эффектов памяти (временная нелокальность) и пространственных корреляций (пространственная нелокальность). Учёт эффектов нелокальности в рамках традиционных подходов приводит к появлению в дифференциальных уравнениях интегрального оператора, где ядро интегрального оператора несёт информацию о природе нелокальности. Для решения таких уравнений интегральные операторы представляются в виде ряда дифференциальных операторов с возрастающим показателем порядка дифференцирования и при наличии малого параметра ограничиваются несколькими членами ряда. В отсутствие малого параметра такой подход оказывается
130
Анализ разностной схемы аналога волнового уравнения ...
Рис. 1. Графики решения задачи (14) при x = 0.7 и различных значениях параметров а и в: 1 — a = 2, в = 2; 2 — a = 2, в = 1.7; 3 — a = 1.7, в = 2
Figure 1. The graphics of solution (14) for x = 0.7 and different values of parameters a and в: 1 — a = 2, в = 2; 2 — a = 2, в = 1.7; 3 — a = 1.7, в = 2
Рис. 2. Графики решения задачи (14) при t = 0.4 и различных значениях параметров а и в: 1 — a = 2, в = 2; 2 — a = 1.7, в = 2; 3 — a = 1.5, в = 2; 4 —
a = 1.4, в = 2
Figure 2. The graphics of solution (14) for t = 0.4 and different values of parameters a and в: 1 — a = 2, в = 2; 2 — a = 1.7, в = 2; 3 — a = 1.5, в = 2; 4 — a = 1.4,
в=2
Рис. 3. Графики решения задачи (14) при t = 0.9 и различных значениях параметров а и в: 1 — а = 2, в = 2; 2 — а = 2, в = 1.7; 3 — а = 2, в = 1.5; 4 —
а = 2, в = 1.4
Figure 3. The graphics of solution (14) for t = 0.9 and different values of parameters а and в: 1 — а = 2, в = 2; 2 — а = 2, в = 1.7; 3 — а = 2, в = 1.5; 4 —
а = 2, в = 1.4
131
В. Д. Бейбалаев, А. З. Якубов
непродуктивным и, кроме того, полученные уравнения также не всегда удаётся разрешить.
СПИСОК ЛИТЕРАТУРЫ/ REFERENCES
1. Ю. И. Бабенко, Метод дробного дифференцирования в прикладных задачах теории тепломассообмена, СПб.: Профессионал, 2009. 584 с. [Yu. I. Babenko, Metod drobnogo differentsirovaniya v prikladnykh zadachakh teorii teplomassoobmena [Method of Fractional Differentiation in Applied Problems of the Theory of Heat and Mass Transfer], St. Petersburg, Professional, 2009, 584 pp. (In Russian)]
2. K. S. Miller, B. Ross, An introduction to the fractional calculus and fractional differential equations, New York, Jon Wiley & Sons. Inc., 1993, xiii+366 pp.
3. I. Podlubny, Fractional differential equations. An introduction to fractional derivatives, fractional differential equations, to methods of their solution and some of their applications, Mathematics in Science and Engineering, vol. 198, San Diego, Academic Press, 1999, xxiv+340 pp.
4. A. A. Kilbas, H. M. Srivastava, J. J. Trujillo, Theory and applications of fractional differential equations, North-Holland Mathematics Studies, Amsterdam, Elsevier, 2006, xv+523 pp.
5. В. В. Васильев, Л. А. Симак, Дробное исчисление и аппроксимационные методы в моделировании динамических систем, Киев: НАН Украины, 2008. 256 с. [V. V. Vasilyev, L. O. Simak, Drobnoe ischislenie i approksimatsionnye metody v modelirovanii dinamicheskikh sistem [Fractional calculus and approximation methods for modeling of dynamic systems], Kiev, NAN Ukraine, 2008, 256 pp. (In Russian)]
6. S. Das, Functional Fractional Calculus, Berlin, Heidelberg, Springer-Verlag, 2011, ix+612 pp. doi: 10.1007/978-3-642-20545-3.
7. В. В. Учайкин, Метод дробных производных, Ульяновск: Артишок, 2008. 512 с. [V. V. Uchaikin, Metod drobnykh proizvodnykh [The method of fractional derivatives], Ulyanovsk, Artichoke Publ., 2008, 512 pp. (In Russian)]
8. А. М. Нахушев, Уравнения математической биологии, М.: Высшая школа, 1995. 301 с. [A. M. Nakhushev, Uravneniya matematicheskoy biologii [Equations of mathematical biology], Moscow, Vysshaya Shkola, 1995, 301 pp. (In Russian)]
9. А. М. Нахушев, Дробное исчисление и его применение, М.: Физматлит, 2003. 271 с. [A. M. Nakhushev, Drobnoe ischislenie i ego primenenie [Fractional calculus and its applications], Moscow, Fizmatlit, 2003, 271 pp. (In Russian)]
10. А. М. Нахушев, Элементы дробного исчисления и их применение, Нальчик: КБНЦ РАН, 2003. 299 с. [A. M. Nakhushev, Elementy drobnogo ischisleniya i ikh primenenie [Elements of Fractional Calculus and Their Application], Nalchik, Kabardino-Balkarsk. Nauchn. Tsentr Ross. Akad. Nauk, 2003, 299 pp. (In Russian)]
11. В. Д. Бейбалаев, “Одношаговые методы решения задачи Коши для обыкновенного дифференциального уравнения с производными дробного порядка” // Вестник Дагестанского государственного университета, 2011. №6. С. 67-72. [V. D. Beybalaev, “One-step methods on the solution of the Cauchy problem for ordinary differential equations with fractional order derivatives”, Vestnik Dagestanskogo gosudarstvennogo universiteta, 2011, no. 6, pp. 67-72. (In Russian)].
12. У. Г. Пирумов, Численные методы, М.: Дрофа, 2004. 224 с. [U. G. Pirumov, Chislennye metody [Numerical methods], Moscow, Drofa, 2004, 224 pp. (In Russian)]
13. Р. Р. Нигматуллин, “Дробный интеграл и его физическая интерпретация”// ТМФ, 1992. Т. 90, №3. С. 354-368; R. R. Nigmatullin, “Fractional integral and its physical interpretation”, Theoret. and Math. Phys., 1992, vol. 90, no. 3, pp. 242-251. doi: 10.1007/ BF01036529.
Поступила в редакцию 01/II/2013; в окончательном варианте — 05/XII/2013; принята в печать — 17/I/2014.
132
Analysis of the Difference Scheme of Wave Equation Equivalent ...
MSC: 35R11, 65M06
ANALYSIS OF THE DIFFERENCE SCHEME OF WAVE EQUATION EQUIVALENT WITH FRACTIONAL DIFFERENTIATION OPERATOR
V. D. Beybalaev1, A. Z. Yakubov2
1 Dagestan State University,
43a, M. Gadzhiev st., Makhachkala, 367025, Russian Federation.
2 Daghestan State Institute of National Economy,
5. D. Ataev st., Makhachkala, 367008, Russian Federation.
Analysis of the difference scheme of boundary-value problem for the wave equation analogue is made. Explicit and implicit difference schemes for numerical solution of the first boundary-value problem for the wave equation analogue with Caputo fractional differentiation operator are investigated, and the stability criteria for these difference schemes are proved by the harmonic Fourier method. Estimates for eigenvalues of the operator of transition from one time layer to another are obtained. Computational experiment on the analysis of the given difference scheme has been performed for the example. The graphs of the numerical solution of the boundary-value problem for the wave equation with the operator of fractional differentiation having different values of parameters of fractional differentiation a and в have been built. Change of the period of fluctuations under transition to a fractional derivative is established. On an example it is shown that parameters a and в become managing directors.
Keywords: wave equation analogue, Caputo fractional partial derivative, difference scheme.
Received 01/II/2013;
received in revised form 05/XII/2013;
accepted 17/I/2014.
ISSN: 2310-7081 (online), 1991-8615 (print); doi: http://dx.doi.org/10.14498/vsgtu1209 © 2014 Samara State Technical University.
Citation: V. D. Beybalaev, A. Z. Yakubov, “Analysis of the Difference Scheme of Wave Equation Equivalent with Fractional Differentiation Operator”, Vestn. Samar. Gos. Tekhn. Univ., Ser. Fiz.-Mat. Nauki [J. Samara State Tech. Univ., Ser. Phys. & Math. Sci.], 2014, no. 1 (34), pp. 125-133. doi: 10.14498/vsgtu1209. (In Russian)
Authors Details: Vetlugin D. Beybalaev (Cand. Phys. & Math. Sci.), Associate Professor, Dept. of Applied Mathematics. Amuchi Z. Yakubov (Cand. Phys. & Math. Sci.), Associate Professor, Dept. of Computer Science & Computer Engineering.
E-mail addresses: [email protected] (V.D. Beybalaev, Corresponding author), [email protected] (A.Z. Yakubov)