CC BY
50УДК 51-7:519.6
Л.И. Мороз
ЧИСЛЕННОЕ РЕШЕНИЕ ОДНОГО КЛАССА НАЧАЛЬНО-ГРАНИЧНЫХ ЗАДАЧ ДЛЯ УРАВНЕНИЯ ДИФФУЗИИ ДРОБНОГО ПОРЯДКА
В работе рассматривается задача для дробно-дифференциального уравнения с частными производными с начальными и граничными условиями, описывающая процесс диффузии. С помощью неявной конечно-разност-ной схемы был получен разностный аналог рассматриваемой задачи. Разработана прикладная программа для численного моделирования процесса аномальной диффузии. Приведены иллюстрации результатов вычислительных экспериментов при варьировании шагов по времени и пространственной координате, проведен анализ полученных результатов.
Ключевые слова: аномальная диффузия, дифференциальное уравнение дробного порядка, аппроксимация, формула Грюнвальда - Летникова, явная схема.
NUMERICAL SOLVING ONE CLASS OF INITIAL-BOUNDARY PROBLEMS FOR FRACTIONAL DIFFUSION EQUATION
The paper considers the initial-boundary problem for the fractional partial differential equation describing the diffusion process. The difference analogue of the problem was described by means of implicit scheme. The program application was developed to numerical simulate anomalous diffusion. Computational experiments results with varying time and spatial coordinate steps are presents. The analysis of the simulation results is also carried out
Key words: anomalous diffusion, fractional-order differential equation, approximation, Grtinwald - Letnikov formula, implicit scheme.
Введение
В настоящее время при исследовании реальных процессов и явлений методами математического моделирования использование строго детерминированного подхода связано с существенными ограничениями. Во многих случаях рассматриваемые физические системы обладают свойствами фрактальности, обусловленными сложным геометрическим строением поверхностей, неоднородностью динамических характеристик и присутствием эффекта наследственности. Теория фракталов нашла применение в описании геометрических свойств сложных объектов, в анализе и прогнозировании поведения динамических систем и процессов. При этом для моделирования динамики процессов и явлений во фрактальных системах часто прибегают к аппарату дробно-дифференциального исчисления. Дробную производную по времени используют, чтобы указать, что рассматриваемый процесс обладает памятью, дробная производная по координате показывает, что процесс протекает в самоподобной неоднородной среде.
Процессы нагрева твердых тел в неравновесных условиях или диффузии примесей в грунте, распространения тепла в высокопористых средах также обладают свойствами, присущими фракталам. Например, в работе [1] описан дробно-дифференциальный подход для моделирования обширного клас-
са задач теории тепломассообмена: теплообмена в инертных средах, изотермических процессов массо-обмена, массообмена на границе контакта движущихся сред, теплообмена при наличии нелинейного тепловыделения. В работе [2] представлена математическая модель аномальной диффузии в форме уравнения в частных производных смешанного порядка. В работе [3] описана модель систем с памятью, обеспечивающая математическую формализацию реальных элементов электрических цепей. Решение задачи о квазистатическом и динамическом растяжении тонкого стержня в модели вязкоупруго-го тела максвелловского типа с дробными производными приведено в работе [4]. В работе [5] получены обобщенные уравнения диффузии и дрейфа носителей заряда в неупорядоченных полупроводниках для дисперсионного переноса. В работе [6] показано, что сегнетоэлектрики являются самоподобными объектами с фрактальной доменной структурой, а результаты исследования фрактальных свойств и вычисления импульса тока переключения поляризации представлены в работах [7, 8].
К сожалению, не все дифференциальные уравнения с производными дробного порядка можно решить с помощью аналитических методов, так как последние разработаны для ограниченного класса задач. В этом случае используют численные методы решения [9].
На сегодняшний день существует несколько определений производной нецелого порядка [10], в числе которых определения:
1) Римана - Лиувилля:
Для численного решения прикладных задач часто используют формулу Грюнвальда -Летникова, позволяющую вводить конечно-разностные аппроксимации дифференциальных задач [9, 10]. Выбирая численный метод для решения задачи, в первую очередь ориентируются на определение, по которому будут аппроксимировать производную, так как от этого напрямую зависит точность вычислений.
Математические модели прикладных задач в физике, химии, экономике, гидрологии и других науках часто сводятся к дифференциальным уравнениям с нецелой производной с начальными и граничными условиями. Аппарат дробно-дифференциального исчисления, используемый для модификации таких математических моделей, позволяет учесть наследственные свойства и фрактальность геометрического строения конкретных объектов. Одной из актуальных задач дробного исчисления является уравнения аномальной диффузии [11].
Многие процессы диффузии вещества, теплопроводности и зарядки описываются подобными математическими постановками. В частности, применительно к области сегнетоэлектрических явлений к дифференциальному уравнению дробного порядка приводит математическая модель переключения поляризации.
Цель настоящей работы - разработка и верификация программного решения для численного моделирования процесса диффузии в краевой задаче для дробно-дифференциального уравнения с частными производными.
Будем рассматривать класс начально-граничных задач для одномерного уравнения диффузии:
(1)
2) Капуто:
(2)
Математическая постановка задачи и вычислительная схема
дг и(х^) , ,даи(х^) , ч
--—- = а(х)--—/- + а(х,п,
ж V / •и '
с начальным условием
и(х^ = 0) = ^/(х) (4)
и граничными условиями первого рода
| и(х = Ь^) = ф2(^)9
где 1<а<2; 0<у<1; 1<х<1; ¿>0; а(х)>0 - коэффициент диффузии; д(х^) - функция источника.
Если зафиксировать порядок производной по координате а, а порядок производной по времени у будет лежать в интервале (ОД) , то уравнение (3) опишет процесс субдиффузии (замедленное блуждание); классическую диффузию наблюдаем при у = 1. Если же зафиксировать у, а 1 < а < 2, то (3) - уравнение супер диффузии (ускоренное блуждание); а = 1 дает классический перенос, а- 2 — уравнение классической диффузии.
Определение производной по Грюнвальду - Летникову имеет вид:
Рг]
Ла ^ и у=0
, /('-Я (б)
У
где
и
а! Г(аг + 1) г-а 0 м ч
- ' а ек; п =-; аек; [х] - целая часть от х; Г(ос] -
Л(а-у)1 Г(у + 1)Г(а-у + 1) А
гамма-функция Эйлера.
Нами используется подход, предложенный авторами работы [12]. Сконструирована неявная конечно-разностная схема, основанная на применении формулы Грюнвальда - Летникова (6).
Данная схема является абсолютно устойчивой и аппроксимирует задачу (3)-(5) с порядком точности
0(Ах + А/)
Для дифференциального уравнения (3) сформулируем его конечно-разностный аналог:
1 -741 п 7+1
Е^Х*"1 (7)
1 1-1 / 1+1 X (М)У п=о г,п ' (Ах)
г
т=0
Рис. 1. Конечно-разностный шаблон _ Т(п-у) _ Г (ти-а)
для неявной схемы. ГДС ~ г(-у)Г(« + 1) ; ^ " Г(-а)Г(/и + 1) " Н°РМИ"
рованные веса Грюнвальда - Летникова.
В нашем случае потенциальный интерес в численной реализации модели с математической постановкой (3)-(5) обусловлен решением прикладной задачи расчета электрических характеристик сегнетоэлектрика при фазовом переходе. Используя термодинамический подход, можно ввести фрактальный аналог уравнения Ландау - Гинзбурга - Халатникова [13]:
= Р + /(А,В,Р) + Е, (8)
где Р - величина поляризации; Е - величина напряженности приложенного поля; 8, А, В - термодинамические параметры.
Задача (8) замыкается заданием начального условия и граничных условий I рода.
Тест-пример
Проведем демонстрацию решения задачи (3)-(5) на следующем тест-примере. Рассмотрим начально-граничную задачу для уравнения диффузии
= + -Л),
а/0
ах1'
Г(1,5)'
(9)
где 0<х<1; ¿>0, с начальным условием г/(х,0) = 0 и граничными условиями первого рода
= 0; =
Аналитическим решением задачи (9) является функция г/(х,^) = х2^.
Для решения класса задач в математической постановке (3)-(5) в ППП МаЙаЬ разработано программное приложение. Входные параметры: шаг по координате Ах и времени концы отрезков х, t; значения порядков производных а, у; коэффициент диффузии а(х); функция источника q(x,t)\ граничные и начальные условия. Выходные параметры: массивы х, матрица
решений и, полученная численным и аналитическим методами; абсолютная погрешность Я; равномерная норма £.
Результаты точного и приближенного решений тест-задачи (8) с шагами по координате и времени Ах = = 0,05 при 0 < г < 1 приведены на рис. 2 а и 2 б.
А б
Рис. 2. Аналитическое решение задачи (9) - а, численное решение этой же задачи, полученное с помощью неявной схемы, - б.
На рис. 3 приведены графики невязок численного решения (9):
К =| и(х^ = 1) - = 1) |, (10)
где и(х^ -1); й(х^ = 1) - точное и приближенное значения функции соответственно.
Абсолютные погрешности Ях вычислены при Ах = /Ы = 0,05 , - при Ах = А/1 = 0,02 , К3 -при Ах = А/ = 0,01 и Я4 - при Ax = At & 0,0067 .
Численное решение задачи (9), полученное с использованием неявной схемы, аппроксимирует аналитическое решение с порядком точности 0(Ах + А^).
С уменьшением шага по времени и координате погрешность приближенного решения стремительно уменьшается. Так, равномерная норма ^ = шах | Я | при при Ах = А/^ = 0,05 составляет 0,004,
1
при уменьшении шагов в 5 раз - Е, « 0,0008308 .
Заключение
В ходе работы была разработана прикладная программа в пакете прикладных программ МайаЬ, предназначенная для численного моделирования диффузионного процесса в постановке краевой задачи для дробно-дифференциального уравнения с частными производными. После анализа полученных результатов можно сделать вывод, что реализованная неявная схема дает приемлемую аппроксимацию и в дальнейшем может быть применена для решения прикладной задачи о моделировании поляризационных характеристик сегнетоэлектрика в процессе переключения поляризации.
1. Бабенко, Ю.И. Метод дробного дифференцирования в прикладных задачах теории тепломассобмена. -СПб.: НПО «Профессионал», 2009. - 584 с.
2. Васильев, В.В. Дробное исчисление и аппроксимационные методы в моделировании динамических систем / В.В. Васильев, JI.A. Симак. - Киев: HAH Украины, 2008. - 256 с.
3. Petras, I. Fractional-order nonlinear systems:modeling, analysis and simulation // HEP, Springer. Nonlinear Physical Science. - 2011. - 235 p.
4. Корчагина, A.H. Использование производных дробного порядка для решения задач механики сплошных сред // Известия Алтайского государственного университета. - 2014. - Т. 1, №1 (81). - С. 65-67.
5. Сибатов, Р.Т. Дробно-дифференциальная кинетика переноса заряда в неупорядоченных полупроводниках / Р.Т. Сибатов, В.В. Учайкин // Физика и техника полупроводников. - 2007. - Т. 41, вып. 3. - С. 346-351.
6. Галиярова, Н.М. Диэлектрическая спектроскопия сегнетоэлектриков, фрактальность и механизмы движения доменных и межфазных границ: Дис. ... д-ра физ.-мат. наук. - Воронеж, 2006. - 399 с.
7. Масловская, А.Г., Барабаш, Т.К. Исследование фрактальных закономерностей процессов переключения поляризации сегнетоэлектрических кристаллов в инжекционном режиме // Поверхность. Рентгеновские, син-хротронные и нейтронные исследования. - 2012. - № 1. — С. 1-8.
8. Мороз, Л.И., Масловская, А.Г. Метод прогноза и коррекции в задаче численного моделирования фрактальной динамики доменных границ сегнетоэлектриков // Вестник АмГУ, серия «Естественные и экономические науки».-2018.-Вып. 83.-С. 11-18.
9. Diethelm, К., Ford, N.J., Freed, A.D., Luchko, Yu. Algorithms for the fractional calculus: a selection of numerical methods // Comput. Methods Appl. Mech. Eng. - 2005. - V. 194. - P. 743-773.
10. Самко, С.Г., Килбас, A.A., Маричев, О.И. Интегралы и производные дробного порядка и некоторые их приложения - Минск: Наука и техника, 1987.-688с.
11. Сластушенский, Ю.В. Математическое моделирование аномальной диффузии с использованием дробно-дифференциальных уравнений и дискретно-элементных моделей: Дис. .. .канд. физ.-мат. наук. - М., 2013. - 105 с.
12. Петухов, A.A., Ревизников, Д.Л. Алгоритмы численного решения дробно-дифференциальных уравнений // Вестник МАИ. - 2009. - Т. 16. - С. 228-234.
13. Физика сегнетоэлектриков: современный взгляд / под ред. K.M. Рабе, Ч.Г. Анна, Ж.-М. Трискона; пер. с англ. - М.: БИНОМ. Лаборатория знаний, 2011. - 440 с.
