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

Вестн. Сам. гос. техн. ун-та. Сер. Физ.-мат. науки. 2012. № 1 (26). С. 157—165

Математическое моделирование

УДК 519.248: 681.5.001.3

ПАРАМЕТРИЧЕСКАЯ ИДЕНТИФИКАЦИЯ ЗАДАЧИ ТИПА КОШИ ДЛЯ ОДНОГО ДРОБНОГО ДИФФЕРЕНЦИАЛЬНОГО УРАВНЕНИЯ

А. С. Овсиенко

Самарский государственный технический университет, 443100, Россия, Самара, ул. Молодогвардейская, 244.

E-mail: [email protected]

Предложен метод параметрической идентификации задачи типа Коши для дифференциальных уравнений, содержащих дробный дифференциальный оператор Римана—Лиувилля порядка а £ (0,1) по мгновенным значениям результатов наблюдений. В основе метода лежит вычисление среднеквадратичных оценок коэффициентов линейно-параметрической дискретной модели функции, аппроксимирующей аналитическое решение. Проведены численно-аналитические исследования, результаты которых позволяют сделать вывод о высокой эффективности предложенного метода.

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

В последние годы усиливается интерес к динамическим системам, описываемым дифференциальными уравнениями, содержащими дробные дифференциальные операторы. Этот факт легко объяснить многочисленными приложениями дробного исчисления в задачах реологии, вязкоупру гости, аномальной диффузии в пористых (фрактальных) структурах, теории автоматического управления, физической химии, биологии и т.д. [1]. В связи с этим встаёт вопрос о параметрической идентификации систем, описываемых дифференциальными уравнениями с дробными производными. Отсутствие известных методов параметрической идентификации таких систем приводит к необходимости разработки эффективных методов, позволяющих определять параметры исходной системы.

Рассмотрим дифференциальное уравнение следующего вида:

П&у -щ = 0, (1)

где

1 d ГТ y(x)dx

мУ Г (1 - a) dt J0 (t - х)а

Анна Сергеевна Овсиенко, аспирант, каф. прикладной математики и информатики.

— дифференциальный оператор Римана—Лиувилля порядка а € (0,1), £ > 0. Для уравнения (1) рассмотрим задачу типа Коши:

11т^-12/ = а1. (2)

Решение задачи (1), (2) можно представить [2] в виде

= а\ Ехр(а, а; /х; ¿), (3)

где Ехр(<т, а] /х; ¿) = Е0-(^сг, а) — двухпараметрическая дробная экспоненциальная функция [3],

п=о

— функция типа Миттаг—Леффлера [4]. Таким образом,

у{1) =а1Г"1Е«(/хГ;а). (4)

Известно [4], что функция типа Миттаг—Леффлера, входящая в (4), при помощи асимптотических формул может быть представлена в виде

Еа(г; а) = -г{1~а)/а ехр (г1/а) + Н(г), Ы -»• оо, (5)

а

причём | а.г§ ^| ^ 7, тта/2 < 7 < тш{7г, 7га},

м _к

ВД-Е^+ОНГ1-")-

С учётом (5) выражение (4) принимает вид

М , ,аЛ-к-1

¡it 01 ^ ^ 7Г

+ о(|/х|-2-м|Га-"м-1). (6)

Выражение (6) содержит параметры /л, сц, а, которые входят в дробное дифференциальное уравнение (1) и его решение (3), и параметр М, определяющий порядок аппроксимирующей функции (6) и погрешность аппроксимации.

Поставим следующую задачу идентификации: при фиксированном значении М и известном значении а вычислить оценки параметров ¡л и сц. Обозначим через ym(t) —функцию, являющуюся аппроксимацией (3) и имеющую вид

ym{t) = ^ ( 1-{1лГ)1/а ехр ((Ю1А*) + £ (/ХГГ"Г(^ + 1)81П(7ГМ) • дt\a 7г /

Отсюда при М = 1 получаем

7Г

(7)

~2Т(а + 1) 8т(7га)

Ыау "3Г(2а + 1) 8т(27га)

Аппроксимация (7) построена на основе асимптотического приближения (6) и является аппроксимацией функции (3), описывающей решение дробного дифференциального уравнения (1) с абсолютной погрешностью порядка 0(|^-3Г2а"1|).

Аппроксимация вида (7) может быть использована при решении задачи параметрической идентификации уравнения (1) в случае, если функция (6) является равномерным асимптотическим разложением решения дифференциального уравнения (1), т.е. отношение двух последовательных слагаемых, входящих в сумму в формуле (6), должно составлять менее одного порядка, что имеет место при выполнении соотношения

/10Г(2а + 1) вт 2тгсЛ 1/а

<0,1, т.е. \ ; .-

\ ¡11 (а + 1) 81П7го; /

Очевидно, что погрешность аппроксимации функции (6) не превышает некоторой величины 5 при £ > / 6. Числен-

но-аналитические исследования предложенной аппроксимации позволяют сделать вывод о том, что выражение (5) является эффективной аппроксимацией функции типа Миттаг—Леффлера начиная с некоторого конечного значения параметра г. На рис. 1 представлены графики точного решения (4) и его аппроксимации, построенной при помощи асимптотической формулы, при значении параметра /л = 1. В этом случае погрешность аппроксимации при £ > ¿о = 1 не превышает 2,5%.

Переобозначая коэффициенты, входящие в выражение (7)

ут{Ь) = с0ерг +

Рис. 1. График точного решения (сплошная линия) и его аппроксимации (штриховая линия)

получаем

где

„ _ . 1/а—1

Со — —V )

а

р =

с 1

а1Г(а + 1) 8ш(7га)

в = —1 — а.

к/л*

(9)

Решим задачу параметрической идентификации процесса, описываемого уравнением (8), при помощи построения линейно-параметрических дискретных моделей (ЛПДМ). Для этого воспользуемся методом, аналогичным описанному в [5]. Полагая в выражении (8) £ = т{к + I), к = 1, 2,... , Ж, где N — объём выборки; I : ¿о = т1, £ > ¿о ф 0; т — период равномерной дискретизации функции ут(Ь), являющейся аппроксимацией точного решения (4) задачи (1), (2). В результате получим дискретный аналог функции ут(Ь). Рассматривая

два последовательных отсчёта дискретной функции у™ = ут(т(к + 1)), после преобразований получаем рекуррентную формулу

уТ = XwT-1 + Ык + 1У + Аз (к + I - 1Г, к = 2,3,..., N.

В виде ЛПДМ и с учётом начальных значений эта формула принимает следующий вид:

{ W = AiC-i + Нк + 1У + A3(fc + I - I)*, к = 2,3,... ,N-, (10)

где Ai = ехр(рт), Л2 = c\ts, Л3 = -AiA2, A4 = А}+гс0 + A2(l + l)s. При использовании среднеквадратичного критерия аппроксимации [6] функции (8) моделью (10) на конечном множестве точек tk = (к + 1)т, к = 1,2,..., N, минимизируется функционал ||е||2 = \\у — у\\2 —> min, или в развернутой форме:

N N

~ Ук)2 = £(ßfc)2 ->■ min, к=1 к=1

где ук — значения функции (3), вычисленные при помощи полученных оценок коэффициентов, у^ — значения, полученные в результате эксперимента.

Вводя в систему (10) случайную аддитивную помеху е в результатах наблюдений, получаем ЛПДМ в виде стохастических разностных уравнений. Будем считать без потери общности, что случайная помеха е имеет нормальное распределение с нулевым математическим ожиданием и некоррелированными значениями в различные моменты времени. С учётом этого ЛПДМ в форме стохастических разностных уравнений можно представить в виде

у? = \А+ец

уТ = ^1Ук-1 + Х2(к + 1У + Нк + 1-1У + г]к, к = 2,3,..., N; (11) Г]к = £к- Ai£fc-1, к = 2,3,... ,N.

В матричной форме ЛПДМ будет иметь вид

Г b = FX + г], \ 7] = Р\£]

где г] = (r]i, г]2, ■ ■ ■, t]n)T — вектор эквивалентного возмущения в стохастическом разностном уравнении, b = (у\, у2,..., 2/w)T ~ вектор правой части, £к = У™ — У™ — случайная помеха в результатах наблюдений, F — матрица регрессоров:

F =

0

У1 У2

о

(2 + I)8

(3 + iy

о

(2 + iy

yN. 1 (N + iy (N + l-iy

а матрица Р\ имеет следующий вид:

Р:

л —

1 0

-А1 1

0 -А

0

О :

0 :

1 :

О -А1

О О

0

1

Оценки коэффициентов А^ (11) могут быть вычислены при помощи итерационной процедуры

д(0) =

№ < = 1.2,...,

сходящейся к решению Л матричного уравнения

рт(рхр1г1рх = рт(рхрУ )~1ъ.

Здесь = — квадратная симметричная матрица разме-

ра N х N.

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

УТ = Л1Со + А2,

УТ = А1УГ-1 + Нк + IУ + Аз (к + I - I)5, к = 2, 3 ... , Ж; У? = Аз,

УТ = А1УГ-1 + л2(к + I)3 - А:А2(к + I - I)5, к = 2,3 ... Ж,

а матрицы регрессоров Р, соответственно, будут иметь вид

^ =

(12) (13)

с0Х\ (1 + 1У 0

У1 (2 +1)3 (1 + 0*

У2 (3 + 1У (2 +0е

ум. 1 (М + 1У (N + 1-1)'

О

У1-А2(1 + 06 У2-А2(2 + 06

О 1

(2 +1)3 О (3 + 1У О

ум-1-\2(М + 1-1у (М + 1У

,-1

О

где Со = — С17г(р5(1 + 8)2Г(—1 — з)8т7г(1 + 8)) и Л2 —оценки, вычисленные на предыдущем шаге по описанной методике.

С помощью итерационной процедуры (13) и с учётом (9) могут быть вычислены оценки параметров модели (8), которые известным образом связаны с коэффициентами Х\, Х2, Аз:

1, „ л- 5 „ Аз — Лг(1 + О5 (л А\ p = z\nXl, С1 = Х2т, с0 =-—щ-. (14)

А1

т

Используя соотношения (9) и (14), вычисляем оценки коэффициентов дробного уравнения (1):

¡1 = ра, а\ = асора~1. (15)

В связи с тем, что количество параметров модели больше, чем число оцениваемых параметров задачи типа Коши (1), (2), проведены дополнительные исследования предложенной процедуры оценивания на основе метода простых итераций. Суть метода заключается в следующем: после вычисления оценок параметров модели при помощи описанного выше метода и пересчёта параметров по формулам (15) вычисляется новое значение параметра С\ на основании выражения (9) и формируется новая итерационная процедура, в которой параметр с\ фиксируется и считается известным, а в процедуре оценивания участвуют только параметры со и р. В этом случае ЛПДМ принимает вид

уТ = х 2,

УТ = ^УТ- V к = 2,3,..., М, (1Ь>

где А1 = ехр(рт), А2 = Соехр(рт(1 + I)). Исследования сходимости данной процедуры показали, что быстрая сходимость позволяет остановить итерационный процесс уже после первого шага. Результаты согласуются с полученными при помощи применения итерационных процедур (12) и (13). Таким образом, использование процедуры (16) позволяет сделать вывод о достоверности полученных результатов.

Проведены численно-аналитические исследования зависимости погрешности оценивания параметров модели (4), а также степени достоверности ЛПДМ-решения (4) от дробного порядка а дифференциального уравнения (1) и периода дискретизации т.

Компьютерное моделирование динамического процесса (4) проводилось при а = 0,5, [I = 1, а\ = 0,7, N = 100, т = 0,02, £ € [1,3] и его результаты представлены на рис. 2-7, на которых приведены результаты оценивания параметров при помощи первой (■) и второй (ж) итерационных процедур. Отклонение модели от точного решения на исследуемом отрезке не превосходит 2,8 %.

На рис. 2, 3 представлены результаты численно-аналитических исследований зависимости погрешности оценивания параметра а\ и среднеквадратичного отклонения в модели от величины случайной аддитивной помехи е. Смещение оценки параметра ц, при 0 < е < 5% не превосходит 0,1%. На рис. 4, 5 представлены результаты численно-аналитических исследований зависимости дисперсий Б оценок коэффициентов ¡лм а\ при наличии случайной аддитивной помехи е в результатах наблюдений. На рис. 6, 7 иллюстрируются результаты численных экспериментов исследования дисперсий Б оценок

коэффициентов ¡л и а\ при наличии случайной аддитивной помехи е в результатах наблюдений.

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

Таким образом, предложен эффективный способ параметрической идентификации дробных дифференциальных операторов порядка а € (0,1), осно-

Ааи, %

А- А"

А % .

О 1 2 3 4 е, %

Рис. 2. Зависимость смещения оценки коэффициента <ц от величины случайной помехи е

5, %

30

20

10

0

¿1

I

АААА А.

ААДА

0 1 2 3 4 е, %

Рис. 3. Зависимость среднеквадратичного отклонения модели от величины случайной помехи е

Оц, %

к А 1

■ ■ ......^ ^ 1 ^...... 1

ч к" ■

11 1 Iй 1

0 1 2 3 4 е, %

Рис. 4. Зависимость дисперсии оценки коэффициента /и от величины случайной помехи е

Эа

30

20

10

0

1 /0 ■ —■

Г4 ■

■ к к"

V ■ Н- "А" А1 ■■ 1...... 1

Рис. 5. Зависимость дисперсии оценки коэффициента <ц от величины случайной помехи е

Dai, %

О

0,2 0,4 0,6 0,8

а

Рис. 6. Зависимость дисперсии оценки коэффициента <ц от порядка дробного дифференциального оператора а

Da?. %

40

30

20

10

0

1 ▲

■ ■ -i ■ А

" "А " ■ ■ iu

......J ■ л

1

0 0,01 0,02 0,03

Рис. 7. Зависимость дисперсии оценки коэффициента сц от периода дискретизации т

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

БИБЛИОГРАФИЧЕСКИЙ СПИСОК

1. Kilbas A. A., Srivastava Н. М., Trujillo J. J. Theory and Applications of Fractional Differential Equations / North-Holland Mathematics Studies. Vol. 204 / ed. J. van Mill. Amsterdam: Elsevier, 2006. 523 pp.

2. Огородников E. H. Некоторые аспекты теории начальных задач для дифференциальных уравнений с производными Римана-Лиувилля // Вестн. Сам. гос. техн. ун-та. Сер. Физ.-мат. науки, 2010. №5(21). С. 10-23. [Ogorodnikov Е. N. Some Aspects of Initial Value Problems Theory for Differential Equations with Riemann-Liouville Derivatives // Vestn. Samar. Cos. Tekhn. Univ. Ser. Fiz.-Mat. Nauki, 2010. no. 5(21). Pp. 10-23].

3. Огородников E.H., Яшагин H. С. Некоторые специальные функции в решении задачи Коши для одного дробного осцилляционного уравнения // Вестн. Сам. гос. техн. унта. Сер. Физ.-мат. науки, 2009. №1(18). С. 276-279. [Ogorodnikov Е. N., Yashagin N. S. Some Special Functions in the Solution To Cauchy Problem for a Fractional Oscillating Equation// Vestn. Samar. Cos. Tekhn. Univ. Ser. Fiz.-Mat. Nauki, 2009. no. 1(18). Pp. 276279].

4. Джрбашян M. M. Интегральные преобразования и представления функций в комплексной области. М.: Наука, 1966. 672 с. [Dzhrbashyan М. М. Integral transforms and representation of functions in the complex domain. Moscow: Nauka, 1966. 672 pp.]

5. Овсиенко А. С., Зотеев В. E. Параметрическая идентификация дробных осцилляторов на основе разностных уравнений / В сб.: Труды шестой Всероссийской научной конференции с международным участием (1-4 июня 2009 г.). Часть 4: Информационные

технологии в математическом моделировании / Матем. моделирование и краев, задачи. Самара: СамГТУ, 2009. С. 61-69. [Zoteev V. Е., Ovsienko A. S. Parametric identification of the fractional oscillator based on difference equation / In: Proceedings of the Sixth All-Russian Scientific Conference with international participation (1-4 June 2009). Part 4 / Matem. Mod. Kraev. Zadachi. Samara: SamGTU, 2009. Pp. 61-69]. 6. Зотеев В. E. Параметрическая идентификация диссипативных механических систем на основе разностных уравнений / ред. В. П. Радченко. М.: Машиностроение-1, 2009. 344 с. [Zoteev V. Е. Parametric identification of dissipative mechanical systems based on difference equations/ ed. V. P. Radchenko. Moscow: Mashinostroenie-1, 2009. 344 pp.]

Поступила в редакцию 13/IX/2011; в окончательном варианте — 20/11/2012.

MSC: 65Р40; 37М05, 26АЗЗ

PARAMETRIC IDENTIFICATION OF CAUCHY PROBLEM FOR ONE FRACTIONAL DIFFERENTIAL EQUATION

A. S. Ovsienko

Samara State Technical University,

244, Molodogvardeyskaya St., Samara, 443100, Russia.

E-mail: [email protected]

The method for parametric identification of Cauchy problem for a fractional differential equation with fractional differential operator of a £ (0,1) degree according to instantaneous values of experimental observations is suggested. The method is based on computation of mean-square estimations for coefficients of linear parametric discrete model of approximation function. Numerically-analytical investigations have been done, the results let us conclude about high efficiency of the method.

Key words: fractional differential operators, parametric identification, linear parametric discrete model, difference equation.

Original article submitted 13/IX/2011; revision submitted 20/11/2012.

Anna S. Ovsienko, Postgraduate Student, Dept. of Applied Mathematics & Computer Science.