Научная статья на тему 'Численное решение обыкновенного дифференциального уравнения в случае дробной разности матричного оператора'

Численное решение обыкновенного дифференциального уравнения в случае дробной разности матричного оператора Текст научной статьи по специальности «Математика»

CC BY
161
26
i Надоели баннеры? Вы всегда можете отключить рекламу.
Ключевые слова
ЧИСЛЕННОЕ МОДЕЛИРОВАНИЕ / ДРОБНЫЙ РАЗНОСТНЫЙ ОПЕРАТОР / ФОРМУЛА ГРЮНВАЛЬДА ЛЕТНИКОВА / COMPUTATION / FRACTIONAL DIFFERENTIAL CALCULUS / GRUNVALD LETNIKOV OPERATOR

Аннотация научной статьи по математике, автор научной работы — Ильин Игорь Александрович, Нощенко Дмитрий Сергеевич, Пережогин Андрей Сергеевич

Рассмотрено численное моделирование решения обыкновенного дифференциального уравнения в случае обобщения разностного оператора. Численные расчеты приведены для уравнения первого порядка. Установлена связь дробной степени матричного оператора СЛАУ с операторами Грюнвальда -Летникова.

i Надоели баннеры? Вы всегда можете отключить рекламу.

Похожие темы научных работ по математике , автор научной работы — Ильин Игорь Александрович, Нощенко Дмитрий Сергеевич, Пережогин Андрей Сергеевич

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.
i Надоели баннеры? Вы всегда можете отключить рекламу.

Computation of differential equation with fractional matrix

Fractional power computation of matrix is represented. Numerical calculations for ordinary differential equation and its fraction analog is done. Results are compared with Grunvald Letnikov operator.

Текст научной работы на тему «Численное решение обыкновенного дифференциального уравнения в случае дробной разности матричного оператора»

УДК 517.3

Численное решение обыкновенного дифференциального уравнения в случае дробной разности матричного оператора *

Ильин И.А.1, 2, Нощенко Д.С.2, Пережогин А.С.1, 2

1 Институт космофизических исследований и распространения радиоволн ДВО РАН, 684034, Камчатский край, c. Паратунка, ул. Мирная, 7

2 Камчатский государственный университет имени Витуса Беринга, 683032, г. Петропавловск-Камчатский, ул. Пограничная, 4

E-mail: [email protected]

Рассмотрено численное моделирование решения обыкновенного дифференциального уравнения в случае обобщения разностного оператора. Численные расчеты приведены для уравнения первого порядка. Установлена связь дробной степени матричного оператора СЛАУ с операторами Грюнвальда - Летникова.

Ключевые слова: численное моделирование, дробный разностный оператор, формула Грюнвальда - Летникова

(с) Ильин И.А., Нощенко Д.С., Пережогин А.С., 2012

MSC 34A08

Computation of differential equation with fractional matrix

Ilyin I.A.1, 2, Noshchenko D.S.2, Perezhogin A.S.1, 2

1 Institute of Cosmophysical Researches and Radio Wave Propagation Far-Eastern Branch, Russian Academy of Sciences, 684034, Kamchatskiy Kray, Paratunka, Mirnaya st., 7, Russia

2 Vitus Bering Kamchatka State University, 683032, Petropavlovsk-Kamchatsky, Pogranichnaya st., 4, Russia

E-mail: [email protected]

Fractional power computation of matrix is represented. Numerical calculations for ordinary differential equation and its fraction analog is done. Results are compared with Grunvald -Letnikov operator.

Key words: computation, fractional differential calculus, Grunvald - Letnikov operator

(c) Ilyin I.A., Noshchenko D.S., Perezhogin A.S., 2012

*Работа выполнена при поддержке Министерства образования и науки РФ в рамках программы стратегического развития ФГБОУ ВПО «Камчатский государственный университет имени Витуса Беринга» на 2012-2016 гг.

Введение

Математическое моделирование физических процессов, протекающих в реальных системах с фрактальными свойствами, выполняется с помощью аппарата дробного дифференцирования [4]. В таких задачах применяются операторы с дробным порядком а, которые в случае целого порядка оператора совпадают с классическими дифференциальными операторами.

В настоящей работе предложен подход к построению модели с дробным оператором через замену классической производной на разностное соотношение и дальнейшее вычисление дробной степени матрицы в системе линейных алгебраических уравнений. Численное моделирование выполнено для обыкновенного дифференциального уравнения первого порядка с начальным условием.

Постановка задачи

Рассмотрим задачу Коши для обыкновенного дифференциального уравнения:

О? + у(г ) = 0, у(0) = 1. (1)

Для этого уравнения запишем разностную схему:

Іу Уі+1 - Уі . „ 1

— «----------, і = 0..п, у0 = 1.

аг т

Таким образом, численную схему для уравнения (1) можно записать в матричном виде:

(і о о о .о Уо (і о о о. .о Уо о

-і і о о .о Уі і о о о. .о Уі о

і о -і і о .о У2 о і о о. .о У2 о

т о о -і і . о Уз + о о і о. . о Уз о

.о. о о о .. .і. W .о. о о о. .. .о. W .о.

(2)

Введем обозначение для матрицы в (2), которая соответствует оператору дифференцирования системы (1):

A=

( і о -і і

о

о

(3)

\о О .. 1/

Вычислим дробную степень а для матрицы А. Положим а = 2. Используем метод Ньютона решения нелинейных уравнений. Вычислительная схема матрицы В = \/Л:

Вп+1 = В+Л'В-1, Во = Л. (4)

2

Используя итерационную схему (4) и систему (2), можно получить численную схему для разностного аналога дробной степени исходного дифференциального оператора уравнения (1).

Реализуем итерационный процесс (4) с начальным приближением Во в виде (3). Алгоритм вычисления запрограммирован в системе компьютерной алгебры Maxima [2]. Приведем пример для матрицы размерности 4x4 (5).

B=

/ 1 0 0 0

-.5 1 0 0

-.125 -.5 1 0

- .0625 -.125 -.5 1

(5)

После этого мы получаем систему (2) с основным оператором представленным матрицей В.

/ 1 0 0 0

1 -.5 1 0 0

-.125 -.5 1 0

- .0625 -.125 -.5 1

У1 +

У2

(т* 0 0 ^ /У0^

0 0 0 У1

0 0 0 У2

0 0 0 W

= 0.

(6)

Решение системы с дробной степенью матрицы А имеет вид, отличный от классического решения системы (2). На рис. 1 приведены решения аналогичных систем с размерностью матриц 100x100, шаг г = 0.05.

Рис. 1. Решение систем линейных уравнений: 1 - классическое решение системы 2-го порядка (2); 2 - линия решение системы порядка 3/2 (6)

Решение системы с матрицей В имеет особенность в начальной точке. При этом асимптотическое поведение решения отличается более длительным режимом затухания в отличии от классического решения.

Сравнение с формулой Грюнвальда - Летникова

Обобщение разностного оператора на случай дробного порядка выполняется с помощью формальной замены целого порядка на дробный. В общем случае, разностный оператор п-порядка с шагом т записывается в виде:

8n f (x) = У Г(п + 1)______________

f () L k!r(n - k + 1)

f (x - kT).

(7)

Полагая порядок оператора п дробным, обозначим через V дробный порядок оператора и получим обобщение [5]:

8v f (x)= У r(v + 1)

f (x) k=0 k!r(v - k + 1)

f (x - kT).

(В)

Корень квадратный из разностного оператора первого порядка (З):

( 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 \

-0.5 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0

-0.125 -0.5 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0

-0.0625 -0.125 -0.5 1.0 0.0 0.0 0.0 0.0 0.0 0.0

-.039063 -0.0625 -0.125 -0.5 1.0 0.0 0.0 0.0 0.0 0.0

-.027344 -.039063 -0.0625 -0.125 -0.5 1.0 0.0 0.0 0.0 0.0

-.020508 -.027344 -.039063 -0.0625 -0.125 -0.5 1.0 0.0 0.0 0.0

-.016113 -.020508 -.027344 -.039063 -0.0625 -0.125 -0.5 1.0 0.0 0.0

-.013092 -.016113 -.020508 -.027344 -.039063 -0.0625 -0.125 -0.5 1.0 0.0

^—0.01091 -.013092 -.016113 -.020508 -.027344 -.039063 -0.0625 -0.125 -0.5 LO.y

1

Hv f (x) = lim — У v

r(v + 1)

т^0 Tv k=0 k!r(v - k + 1)

f (x - kт)

Проведем вычисления в пакете Maxima для коэффициентов входящих в формулу (В) с показателем дробности v = 2

nu : 1/2;

eq : (-1^k*gamma(nu+1)/gamma(k+1)/gamma(nu-k+1);

makelist( subst( j, k, eq ),j,0,10),eval,float;

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

В результате вычисления для коэффициентов дробного разностного оператора получаем значения, которые совпадают со значениями матрицы B:

[1, -.5, -.125, -.0625, -.039063]. (9)

Получаем согласование вычисления квадратного корня из матрицы A c формулой (В), которая получена формальной процедурой обобщения разностного оператора на случай дробного порядка.

Возьмем определение оператора Грюнвальда - Летникова Hv и установим связь с вычисленным корнем квадратным из матрицы A:

(10)

Если не вычислять предел при т ^ 0, то Н2 совпадает с выражением корня квадратного из матрицы А. При этом ряд, входящий в формулу Грюнвальда - Летникова, обрывается при задании начального условия. Полученный результат говорит о том, что вычисление значения в текущей точке зависит от всех предыдущих значений.

Матричный оператор степени больше чем единица

Умножим матрицу (3) на (5) в случае размерности 10x10, и в результате матричный оператор будет иметь вид:

AB =

1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0

-1.5 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0

0.375 -1.5 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0

0.0625 0.375 -1.5 1.0 0.0 0.0 0.0 0.0 0.0 0.0

.02344 0.0625 0.375 -1.5 1.0 0.0 0.0 0.0 0.0 0.0

.01172 .02344 0.0625 0.375 -1.5 1.0 0.0 0.0 0.0 0.0

.006836 .01172 .02344 0.0625 0.375 -1.5 1.0 0.0 0.0 0.0

.004395 .006836 .01172 .02344 0.0625 0.375 -1.5 1.0 0.0 0.0

.003021 .004395 .006836 .01172 .02344 0.0625 0.375 -1.5 1.0 0.0

^.002182 .003021 .004395 .006836 .01172 .02344 0.0625 0.375 -1.5 1.0.

Вычислим коэффициенты входящие в определение оператора Грюнвальда - Лет-никова при а = |. Получим значения коэффициентов:

[1.0, -1.5,0.375,0.0625, .02344, .01172, .006836, .004395, .003021, .002182, .001637].

Коэффициенты в точности совпадают с результатом, который получен в произведении А ■ В.

Рассмотрим задачу Коши для обыкновенного дифференциального уравнения второго порядка:

Л2у

+ y(t )= 0, y(0) = 1, y'(0) = 0.

(11)

В качестве аппроксимации оператора второго порядка выберем вторую разность:

Л 2у Уі+1 - 2у,- + у-1 .

*2 « --------т2-----------------, ' = °"“, (12)

тогда систему линейных уравнений для уравнения (11) можно записать в виде:

1 0 0 0 /У>\ - __L т2 0 0 0 /У>\

-2 1 0 0 У1 + 1 т2 0 0 0 У1

1 -2 1 0 У2 0 1 0 0 У2

0 1 -2 1 У3 V 0 0 1 0 У3

= 0. (13)

В качестве вычислительного эксперимента зададим значение а = | и система

линейных алгебраических уравнений примет вид:

3

т 2

1 0 0 0

-1.5 1 0 0

.375 -1.5 1 0

^.0625 .375 -1.5 1

/W\

У1

У2

+

(-\ 0 0 0\

т 2

1 0 0 0

2т 3/2 0

1 0 0

/V\

У1

У2

= 0.

Построим графики решений для системы второго порядка и - (рис. 2).

1

2

1

Рис. 2. Решение систем линейных уравнений: 1 - решение системы 2-го порядка (13); 2 - решение системы порядка 3/2 (14)

Как и в случае системы для системы с порядком дробности а = - имеется особенность в начальной точке. Однако, асимптотический характер решения совпадает с тем, который хорошо известен для дробного осциллятора [3]. Подобная динамика решения хорошо согласуется с задачами Коши для обыкновенных уравнений с дробной производной Римана - Лиувилля. Чтобы устранить особенность в начальной точке, можно использовать методику предложенную в работе [1]. В целом, представленные результаты могут быть использованы при численном моделировании динамики систем с памятью. В частности, в теории вязко-упругости можно выполнить моделирование релаксационных процессов для сложных природных сред и композитных материалов.

Заключение

Полученные численные решения имеют принципиальное отличие от экспоненциального решения за счет более медленного спадания решения на асимптотике. Характер решения указывает на особенности динамики таких систем, которые проявляются в наличии эффектов памяти. Связь, установленная между разностным оператором дробного порядка Грюнвальда - Летникова и дробной степенью матрицы, указывает на прямую связь двух подходов.

Библиографический список

1. Ильин И.А., Нощенко Д.С., Пережогин А.С. О дробной степени разностного оператора для обыкновенного дифференциального уравнения // Математическое моделирование фрактальных процессов, родственные проблемы анализа и информатики: материалы II Междунар. конф. молод. ученых. Нальчик; Эльбрус, 2012. С. 107-109.

2. Электронный ресурс: URL: maxima.sourceforge.net/

3. Нахушев А.М. Дробное исчисление и его применение. М.: Физматлит, 2003. 272 с.

4. Тарасов В.Е. Модели теоретической физики с интегро-дифференцированием дробного порядка. М., Ижевск: РХД, 2011. 568 с.

5. Учайкин В.В. Метод дробных производных. Ульяновск: Артишок, 2008. 512 с.

Поступила в редакцию / Original article submitted: 03.12.2012

i Надоели баннеры? Вы всегда можете отключить рекламу.