CC BY
23
МАТЕМАТИКА И МЕХАНИКА
УДК 517.927
Исследование численного метода решения краевой задачи для дифференциального уравнения с дробной производной по времени*
Н.Б. Алимбекова1,2, Д.Р. Байгереев2, М.Н. Мадияров2
'Казахский Национальный педагогический университет им. Абая (Алматы, Казахстан)
2Восточно-Казахстанский государственный университет им. С. Аманжолова (Усть-Каменогорск, Казахстан)
Study of a Numerical Method for Solving a Boundary Value Problem for a Differential Equation with a Fractional Time Derivative
N.B. Alimbekova1,2, D.R. Baigereyev2, M.N. Madiyarov2
'Abai Kazakh National Pedagogical University (Almaty, Kazakhstan)
2Sarsen Amanzholov East Kazakhstan State University (Ust-Kamenogorsk, Kazakhstan)
В настоящее время замечается повышенный интерес к проблеме численной реализации моделей многофазной фильтрации в связи с ее огромной экономической значимостью в нефтедобывающей промышленности, гидрологии и управлении ядерных отходов. В отличие от классических моделей фильтрации, модели фильтрации в сильнопористых трещиноватых пластах с фрактальной геометрией скважин изучены недостаточно полно. Решение данной задачи сводится к решению системы дифференциальных уравнений с дробными производными. Построена конечно-разностная схема для решения начально-краевой задачи для уравнения конвекции-диффузии с производной дробного порядка по времени в смысле Капуто-Фабрицио. Получены априорные оценки для решения разностной задачи в предположении существования решения задачи в классе достаточно гладких функций, которые доказывают единственность решения и устойчивость разностной схемы. Показана сходимость решения разностной задачи к решению исходной дифференциальной задачи со вторым порядком по временной и пространственной переменным. Представлены результаты вычислительных экспериментов, подтверждающие достоверность теоретического анализа.
Ключевые слова: уравнение дробного порядка, производная Капуто-Фабрицио, устойчивость, сходимость, априорная оценка.
DOI 10.14258/izvasu(2020)4-10
Recently, there has been an increased interest in the problem of numerical implementation of multiphase filtration models due to its enormous economic importance in the oil industry, hydrology, and nuclear waste management. In contrast to the classical models of filtration, filtration models in highly porous fractured formations with the fractal geometry of wells are not fully understood. The solution to this problem reduces to solving a system of differential equations with fractional derivatives. In the paper, a finite-difference scheme is constructed for solving the initial-boundary value problem for the convection-diffusion equation with a fractional time derivative in the sense of Caputo-Fabrizio. A priori estimates are obtained for solving a difference problem under the assumption that there is a solution to the problem in the class of sufficiently smooth functions that prove the uniqueness of the solution and the stability of the difference scheme. The convergence of the solution of the difference problem to the solution of the original differential problem with the second order in time and space variables is shown. The results of computational experiments confirming the reliability of theoretical analysis are presented.
Key words: fractional differential equation, Caputo-Fab-rizio fractional derivative, stability, convergence, a priori estimate.
*Работа выполнена при поддержке грантового финансирования научно-технических программ и проектов Министерства науки и образования Республики Казахстан, ИРН АР08053189, 2020-2022 годы.
Уравнения дробного порядка находят применение при описании большого класса физико-химических процессов, протекающих в средах с фрактальной геометрией, а также при математическом моделировании экономических и социально-биологических явлений. Одним из важных примеров применения уравнений данного типа являются уравнения, описывающие течение многофазной жидкости в сильнопористых трещиноватых пластах с фрактальной геометрией скважин. Поэтому разработка аналитических и вычислительных методов решения уравнений дробного порядка является актуальной и важной проблемой.
ссследованию вопросов существования и единственности решения задач для уравнений данного типа посвящено множество работ [1-5]. Разработке численных методов решения краевых задач для уравнений с производной дробного порядка посвящены работы [2,3,6-10].
В данном исследовании построен численный метод решения уравнения с дробной производной по времени в смысле Сапуто-Фабрицио. Выведены нормализованные веса в аппроксимационной формуле для дробной производной. Построена неявная конечно-разностная схема второго порядка аппроксимации по времени и по пространственной переменной. Доказана устойчивость предложенной схемы, а также сходимость со скоростью, равной порядку аппроксимации. Полученные результаты подтверждаются численными расчетами, проведенными для тестовой задачи.
1. Постановка задачи. В области Qт = {(х,4) : 0 < х < 1, 0 < 4 < 1} рассмотрим задачу
д
ди д
ди
даи
да
1
1 — а
и' (т) dт
ем—)
^ 1 =
1а
Предполагается, что для коэффициентов и правой части уравнения (1) выполняются условия
k(х,г) е с1-0 (дт), Г (х,1),1 (х,1) е с ^), (4)
с1 < k (х,Ь) < С2, ^ (х,Ь)\ < С2,
С1 > 0, С2 > 0, 2С1 > с2.
(5)
Предположим, что существует решение задачи (1)-(3) в классе достаточно гладких функций.
2. Постановка разностной задачи. Для
численного решения задачи (1)-(3) применим метод конечных разностей. В области [0,1] х [0, Т] введем равномерную конечно-разностную сетку
шьт = {(х,п) : х= Л, Ьп = пт,
i = 0Ж п = 0^, xN = 1, № = Т} .
Иа сетке дифференциальной задачи (1)-(3) поставим в соответствие разностную схему с весами порядка аппроксимации О (Н2 + т2):
Д&„+„ Уi = + пи+й+УаЬ
+ (Ы/Л + Фп, (х,1) е , (6)
\ / х,0
у0 = Pi,
у(0а) = у№г) =0, п> 0,
(7)
(8)
где Д^ ^ у - дискретный аналог дробной производной Сапуто-Фабрицио порядка а, 0 < а < 1:
да _ ' \ ' а,а s
Д0tn+стУ = а gn-syt ;
коэффициенты даа,а определяются по формуле
90 = А0
при п = 0,
Ааа+ва'а,
s = 0,
да = qМ дхх + ^(хдТх) + 1 (хЛ, (1)
0 < х < 1, 0 <г < 1,
и (х, 0) = р (х), (2)
и (0,4) = и (1,4) = 0, (3)
где 0 < а < 1. Производная дробного порядка определена в виде
t
а = { Аа,.+ва+1 - ва>°, 1 < s < п -1,
Ап вп ,
при п > 0, где
1 — е-^та О = 0 ла,а _ I 1 е , о 0,
Аё = ^ ^-1 О> 1,
е~г
-(,<*+') ,
е1Т (1т - 2)+ 1т + 2 1
в, = ---7—;—N-, О > 1,
а также введены обозначения
п- = \(п -\п\), п+ = 22(п + \п\),
£п = М хо- 1,4
г+а
п = q (х,п+а) , П0 k (хо,гп+а) ,
фП = I Сх,п+а), а = 1 - а, 9аа > 0,
у(а) = ауп+1 +(1 - а) уп.
а
3. Устойчивость и сходимость разностной схемы. Докажем несколько вспомогательных лемм.
Лемма 1. Для любой функции у (£), заданной на сетке , справедливо неравенство
У(а)д^гп+,у > у2.
Данная лемма доказывается аналогично лемме 1 из работы [11]. Ниже буквами ц с индексами будут обозначаться положительные числа, не зависящие от h и т.
Лемма 2. При выполнении условий (5) для решения разностной задачи (6)-(8) справедливо неравенство
1д°
2 дог,
\\уг + т
<
< (а) ||уп+1\\2 + М3 (а) \\уп\\2 + М2 . Доказательство. Умножим уравнение (6) ска-
лярно на у
2До
Ы? + ЧЧ
< 14
+
Теорема 1. При выполнении условий (5) существует то такое, что при т < т0 для решения разностной задачи (6)-(8) справедлива оценка
,,п+1 \
I 2 /и о \ \2 ип
ну1 1 <мм 11у ii + та5 \\м
4 0<т<и
(Д^ у= (V еу{°\у(а)) + + (п+^Ч^У*) +
+ ((^)х ,У(а)) + (м,У(а)) . (10)
Скалярное произведение в левой части тождества оценим с помощью леммы 1:
у,у(^)) >1 (1, даи+„у2) = 2да(„+ Ы.
(11)
Оценим остальные слагаемые с помощью разностного аналога первой формулы Грина, неравенства Коши с е и условий (5), (8):
из которой следует единственность и устойчивость разностной схемы (6)-(8) по начальным данным и правой части.
Теорема 2. При выполнении условий теоремы 1 решение разностной задачи (6)-(8) сходится к решению дифференциальной задачи (1)-(3) и выполняется оценка
НуП+1 - иП+1Н < ^2 + т2) .
Доказательство. Для доказательства сходимости разностной схемы (6)-(8) рассмотрим задачу для разности z = у — и:
д04„+„ z = Сг'^х} + гН+Сг+12'Х} +
+ (Cz*) )х + (М) е , (13) z0 =0, (14)
1°) = л?)
¿о = ¿и
0, п> 0, (15)
(де ф? = мП — д&п+„и? + Пп-^ + +
^ПиХ,}^ . Для решения разностной задачи (13)-(15) выполняется оценка
ип+1\\ < щ> тах \\фш\\
0< ш<п
(16)
(12)
Используя определение у(а) в (9), из (12) получим утверждение леммы.
Лемма 3 [2]. Пусть неотрицательные последовательности у? и М?, п = 0, 1, 2,... удовлетворяют неравенству
даМп+„Уп < Х1У?+1 + х2Уп + Мп, п > 1,
где А1 > 0, Х2 > 0 - некоторые константы. Тогда существует такое то, что при т < то выполняется неравенство
Уп+1 < 1М>[ У0 + тах М
0<ш<п
На основании лемм 2 и 3 доказывается следующая теорема.
где \\фп\\ = О (Н2 + т2). Из априорной оценки (16) следует сходимость решения разностной задачи (6)-(8) к решению дифференциальной задачи (1)-(3).
4. Реализация разностной схемы и анализ вычислительных экспериментов. При
реализации разностной схемы (6)-(8) на каждом временном слое требуется решить систему линейных уравнений
АгУп+11 — Спуп+1 + Вуп+ = i = 1,2, ...,М — 1,
где
в = + 4 (пп + К \)) ,
а}а
с = ^ + « №+1 (пп + \ пГ \ ) —
— ^ (пп — \ пп\)) + « ^ + £+0
2
2
2
2
2
2
Рис. 1. График приближенного решения задачи на различных временных слоях
г-1
Fi =- \coyn -У ]еп-з Ы+1 - У
3=0
+
+п - \п?\) (у? - Уи) +
+(пп + \п?\) (уг+1 - у?) +
+ 1- (£+1 (уп+1 - у?) - ^ (у? - уп-1)) + v?.
Данная система решается методом скалярной прогонки. Непосредственной проверкой можно убедиться в том, что условие устойчивости метода прогонки выполняется.
Для проверки точности разностной схемы (6)-(8) проведен ряд вычислительных экспериментов на примере модельной задачи.
Задача. Рассмотрим уравнение
д
ди 2
= + ~ х2(1 - х) ( г - - -
дЬа дх а V 7
1
7е1*)
-г2 (2 - 4х - 3х2), о <х< 1, о <г< 0.1
с однородными начальным и граничными условиями. Точное решение задачи имеет вид
и (х, г) = г2х2 (1 - х).
При анализе зависимости порядка погрешности от пространственного шага значение временного шага выбрано равным т = 10-5. Значение шага по пространственной переменной Н варьировалось между Н = 10-2 и Н = 10-5. Порядок дробной производной выбран равным а = 0.1, а = 0.45 и а = 0.85. Вычисления проводились до достижения временного слоя п = 1000.
Величина погрешности определялась по формуле
Е = тах тах \у? — и (х,А?
При анализе зависимости порядка погрешности от временного шага значение пространственного шага выбрано равным Н = 10-4. Значение временного шага варьировалось между т = 10-5 и т = 10-8. Порядок дробной производной выбран равным а = 0.3, а = 0.45 и а = 0.85.
В таблицах 1 и 2 приведены значения погрешностей для различных значений параметров а, Н и т.
а
Рис. 2. Приближенное решение задачи
Таблица 1 Анализ погрешности для задачи
Таблица 2 Анализ погрешности для задачи
Н а = 0.95 (а = 0.1) а = 0.775 (а = 0.45) а = 0.575 (а = 0.85)
1 100 1.37 • 10-7 1.09 • 10-7 5.93 • 10-8
1 ,500 5.55 • 10-9 4.58 • 10-9 2.63 • 10-9
1 1000 1.43 • 10-9 1.29 • 10-9 8.67 • 10-10
1 2000 4.00 • 10-10 4.85 • 10-10 4.32 • 10-10
1 5000 1.14 • 10-10 2.65 • 10-10 3.12 • 10-10
10000 7.06 • 10-11 2.35 • 10-10 2.96 • 10-10
т а = 0.85 (а = 0.3) а = 0.775 (а = 0.45) а = 0.575 (а = 0.85)
10-5 1.38 • 10-9 1.29 • 10-9 8.67 • 10-10
10-® 1.15 • 10-11 8.65 • 10-12 2.39 • 10-12
10-7 8.69 • 10-14 4.62 • 10-14 4.63 • 10-15
10-8 5.87 • 10-16 2.14 • 10-16 7.50 • 10-18
На рисунках 1 и 2 приведен график приближенного решения задачи на различных временных слоях.
Из результатов, приведенных в таблице 1, следует, что погрешность является величиной порядка О (Н2). Аналогично из результатов, приведенных в таблице 2, следует, что погрешность является величиной порядка О (т2). Таким образом, вычислительные эксперименты подтвердили, что разностная схема сходится со вторым порядком как по пространственной, так и по временной переменным.
Заключение. Таким образом, в работе построена неявная конечно-разностная схема для дробного дифференциального уравнения с переменными коэффициентами, содержащего дробную производную по времени в смысле Капуто-Фабрицио. Установлены устойчивость и оценка погрешности разностной схемы. Эмпирическая сходимость хорошо согласуется с теоретическими оценками. Полученные результаты будут использованы при построении конечно-элементных методов решения задачи для дифференциального уравнения с дробной производной по времени в смысле Капуто-Фабрицио.
Библиографический список
1. Berdyshev A., Cabada A., Turmetov B. On solvability of some boundary value problem for polyharmonic equation with boundary operator of a fractional order // Applied Mathematical Modelling. 2015. T. 4. DOI: 10.1016/j.apm.2015.01.006.
2. Alikhanov A.A. A new difference scheme for the time fractional diffusion equation // Journal of Computational Physics. 2015. T. 280. DOI: 10.1016/j.jcp.2014.09.031.
3. Berdyshev A., Eshmatov B., Kadirkulov B. Boundary value problems for fourth-order mixed type equation with fractional derivative // Electronic Journal of Differential Equations. 2016. № 36.
4. Agarwal P., Berdyshev A., Karimov E. Solvability of a Non-local Problem with Integral Transmitting Condition for Mixed Type Equation with Caputo Fractional Derivative // Results in Mathematics. 2017. DOI: 10.1007/s00025-016-0620-1.
5. Бештоков М.Х. Нелокальные краевые задачи для уравнения соболевского типа с дробной производной и сеточные методы их решения // Математические труды. 2018. T. 21, № 2. DOI: 10.17377/mattrudy.2018.21.203.
6. Beshtokov M. Boundary value problems for degenerate and degenerate fractional order differential equations with non-local linear source
and difference methods for their numerical implementation // Ufimskii Mathematicheskii Zhurnal. 2019. T. 11, № 2.
7. Kanwal A., Phang C., Iqbal U. Numerical Solution of Fractional Diffusion Wave Equation and Fractional Klein-Gordon Equation via Two-Dimensional Genocchi Polynomials with a Ritz-Galerkin Method // Computation. 2018. T. 6, № 40. DOI: 10.3390/computation6030040.
8. Jin B., Lazarov Y., Liu Y., Zhou Z. The Galerkin finite element method for a multi-term time-fractional diffusion equation // Journal of Computational Physics. 2015. T. 281. DOI: 10.1016/j.jcp.2014.10.051.
9. Morales-Delgado V.F., Gomez-Aguilar J.F., Taneco-Hemandez M.A. Analytical solution of the time fractional diffusion equation and fractional convection-diffusion equation, Revista Mexicana de Fisica. 2019. T. 65. DOI: 10.31349/RevMexFis.65.82.
10. Liu F., Zhuang P., Burrage K. Numerical methods and analysis for a class of fractional advection-dispersion models // Computers and Mathematics with Applications. 2012. T. 64 DOI: 10.1016/j.camwa.2012.01.020.
11. Alikhanov A.A. A priori estimates for solutions of boundary value problems for fractional-order equations Differential Equations. 2010. T. 46. DOI: 10.1134/S0012266110050058.
