Вычисление скорости роста вектора состояний для одной модели обобщенной линейной стохастической системы Текст научной статьи по специальности «Математика»

Н. К. Кривулин

ВЫЧИСЛЕНИЕ СКОРОСТИ РОСТА ВЕКТОРА СОСТОЯНИЙ ДЛЯ ОДНОЙ МОДЕЛИ ОБОБЩЕННОЙ ЛИНЕЙНОЙ СТОХАСТИЧЕСКОЙ СИСТЕМЫ*

1. Введение. При исследовании технических, экономических, производственных и других систем находят применение обобщенные линейные динамические модели вида

z(k) = A(k) <g) z(k — 1),

где A(k) —случайная матрица системы, z(k) —вектор состояний системы, ® —операция умножения матриц, заданная в некоторой идемпотентной алгебре [1, 2].

При анализе реальных систем часто представляет интерес определение средней скорости роста вектора состояний системы

А= ,Иш гИ*)Н>

к—k

где || • || —некоторый идемпотентный аналог обычной векторной нормы.

Во многих случаях доказательство существования указанного предела может быть проведено без труда (например, на основе эргодической теоремы из работы [3]). В то же время вычисление значения предела обычно оказывается довольно трудной задачей даже в случае весьма простых моделей систем. К имеющимся результатам в этой области можно отнести решение задачи вычисления предела [4] для системы второго порядка с матрицей, все элементы которой независимы и имеют одинаковое экспоненциальное распределение с единичным средним. В работе [5] рассматривалась система второго порядка с симметричной матрицей, у которой диагональные элементы независимы и имеют экспоненциальное распределение с единичным средним, а все недиагональные элементы равны нулю. В [6] найдено общее решение для случая системы с треугольной матрицей порядка n при условии, что случайные элементы матрицы имеют распределения вероятностей с конечными средним и дисперсией и могут быть зависимыми.

Целью настоящей работы является решение задачи определения средней скорости роста вектора состояний для системы второго порядка с матрицей, элементы которой независимы и экспоненциально распределены. При этом предполагается, что диагональные элементы имеют одинаковые распределения с параметром ^, а недиагональные элементы — распределения с параметром v = ^. Предложенная модель системы является некоторым обобщением моделей, изученных в [4, 5].

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

* Работа выполнена при финансовой поддержке РФФИ (грант №06-01-00763).

© Н. К. Кривулин, 2007

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

2. Линейная стохастическая динамическая система. Рассмотрим динамическую систему, эволюция которой описывается уравнением

а знак <8> обозначает операцию умножения матрицы на вектор в идемпотентной алгебре со скалярными операциями обобщенного сложения ® и умножения <8>, определенными соответственно как операция максимума и сложение (см., например, [1, 2]).

Векторное уравнение (1) можно представить в виде системы скалярных уравнений

Будем предполагать, что каждая из последовательностей {ак}, {вк}, {1к} и } состоит из независимых и одинаково распределенных случайных величин, а также, что ак, вь 7т и 5п независимы при любых к,1,т,и.

Пусть диагональные элементы ак и 5к матрицы А(к) имеют общее экспоненциальное распределение вероятностей с параметром ^, а недиагональные элементы вк и 7к — экспоненциальное распределение с параметром V. Обозначим функции распределения указанных элементов через

а соответствующие им функции плотности через /(і) и д(і).

Ниже будем предполагать, что ц = V, так как при условии ц = V рассматриваемая система легко сводится к системе, для которой ц = V =1 [4, 5].

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

где обозначение нормы понимается в смысле ||г(к)|| = х(к) ® у(к) = шах{х(к), у(к)}.

Нетрудно показать, например, с использованием эргодической теоремы [3], что в случае, когда элементы матрицы А(к) неотрицательны и имеют конечные средние, указанный предел существует с вероятностью 1, а также существует предел

z(k) = А(к) ® г (к — 1),

(1)

где

х(к) = ак ® х(к — 1) © вк <8> у(к — 1), у(к) = 7к ® х(к — 1) © 5к <£> у(к — 1)

или, с использованием обычных обозначений, в виде

х(к) = шах(«к + х(к — 1), вк + у(к — 1)}, у(к) = шах(7к + х(к — 1), 4 + у(к — 1)}.

^(і) = шах{0,1 — е м4}, О(і) = шах{0,1 — е и1},

Ит ^Е||г(*:)|| = А.

к—— ^ к

Для нахождения предела введем так же, как и в работе [4], новые переменные

Кроме того, заметим, что в силу очевидной симметрии представления системы для всех £ имеет место равенство —^(— £) = —&(£).

Предположим, что последовательность функций —^ при к ^ то сходится к некоторой функции плотности — равномерно по £ (ниже будет показано, что такое предположение для рассматриваемой системы является справедливым).

Тогда будет существовать предел

Обозначив функцию плотности, соответствующую функции распределения Ф, через ф, приходим к выражению для вычисления Л в виде

X (к) = х(к) - х(к - 1), У (к) = у (к) - х(к)

и заметим, что х(к) = X(1) + • • • + X(к).

Тогда уравнения системы можно переписать в виде

X(к) = шах{а, вк + У(к - 1)}, У (к) = шах(7к ,5к + У (к - 1)} - шах{«к,вк + У (к - 1)}.

Введем функции распределения

Фк(і) = Р{Х(к) < і}, *к(*) = Р{У(к) < і},

а соответствующие им функции плотности обозначим через фи(^) и —^(£). Из первого уравнения полученной системы следует, что

а также конечные пределы

Ііш Е[Х(к)] = а, Ііш Е[У(к)] = Ь, Ііш Е[шах{У(к),0}] = Ь+.

Рассмотрим величину

к

|| г (к) У = шах{х(к), у (к)} = шах{У (к), 0} + X (і).

Для средней скорости роста вектора состояний системы Л имеем

11 1 .к

А= Ит —Е||г(/г)|| = Ит —Е[тах{У(/г), 0}] + Ит — Е[Х (г)] = а.

і=1

4. Рекуррентное уравнение для функции плотности —к. По формуле полной вероятности для любого к = 2, 3,... выполняется равенство

K(t, s)-fc_i(s)ds, (3)

где

K(t, s) = P{Y(k) < t|Y(k — 1) = s} = P(max(7k, 5k + s} — max{ak, вк + s} < t}.

Рассмотрим функцию K. Учитывая, что K(t, s) при каждом s является функцией распределения разности случайных величин

£(s) = max{Yfc, Sk + s}, n(s) = max{afc, вк + s},

для которых P{£(s) < t} = G(t)F(t — s) и P{n(s) < t} = F(t)G(t — s), имеем

K(t, s) = / G(u + t)F(u + t — s)(f (u)G(u — s) + F(u)g(u — s))du,

J r

где

{—t, если t < 0, s < 0,

s — t, если t < 0, s > 0,

0, если t > 0, s < 0,

s, если t > 0, s > 0.

Найдем функцию K. При t < 0, s < 0 вычисление интеграла дает

JT(i, s) = -i/e*“ f —-----------—е**Л + -e"4 f 1 - 7---------^-ге^Л e"®-

2 YyU + г/ 2yU, + г/ J 2 \ (/x + z/)(/x + 2z/) у

- ^e(ii+,y)i Г—:—e"“ I e'

2 V M + 2v 2^ + v

При t < 0, s > 0 имеем

s) = ^ f1 - 7--------------e~MS + f —----------

2 \ (p + z/)(2/x + г/) у 2 YyU + г/ ц + 2г/ у

_ [ ___?____________i:__I p-мя

2M + ^ М + 2г/

при t > 0, s < 0 имеем

K{t, s) = 1 - -e-*“ f 1 - 7^-------------------------е"Л e*“ - ( —-------------- --е"Л

2 \ (p + v) (2/x + г/) у 2 YyU + г/ ц + 2г/ у

-Me^(M+,y)i Г_____-_________________-_с

2М \2n + v n + 2v

при і > 0, в > 0 имеем

КЦ, е) = 1- ( —---------------------—--е-^Л - (1 - 7---------------------^-—е-"8 ) е-и8 +

2 \/х + V ¿¡л + V у 2 у (р + */)(а* + 2г/)

_І_ Іг/е-СИ-")* ^__?__________і___р-мя

^ + 2^ 2^ + V

Введем следующие функции

А1{ґ)=1-е-^\, А<2{ґ)=1-е-^, А3{ґ)=1-е-^ І(І,

яіСО = “2(і) = ^е-"1*1, аз (і) = ^(м + ^е4^1*1,

а также функции

Ь1(в) =

^ ~ ' ’ ЄСЛИ 5 ^ °’ е~М8 (: - ы+^+,)е^8). если 5 > °;

ъ Ы = , ^ Iі - (М+.)(М+2.)ЄМ8) - ЄСЛИ 5 ^ °-2 М - ¡¡Т2ЇЇЄ_

- ¡¡Т2їїе "8) > если 5 > °;

-іує1'3 (—^- - ^—е^3) , если в < О, V /+2^ 2^+^ ' 1 — 1

3( ) ^ (її^ - ЛТ2їе 1/а) > если 5 > °-

-^в-^8

Теперь функция К может быть представлена в виде N I К.(£, в), если £ < 0, ^

КЛ , „ К1(£,в)^^3 А4(£Мв).

11 — К1(£, -в), если £> 0; 4=1

Уравнение (3) принимает вид

{3 оо

ЕМ^) / Ь»(в)—й-1(в)^в, если £ < 0,

4=1 з -0 о

1 - £ А4(£) / Ь4(-в)—й-1(в)^в, если £ > 0.

4=1 —О

Учитывая, что —й_1(-в) = —&_1(в), имеем равенство

/ОО /*О

&*(-«)—й_1(«)^в = Ь»(в)—Й_1(в)йв,

О _О

из которого следует, что уравнение (3) можно записать так:

*к(і) = 1 І=1

ЕАі(і) / Ьі(5)^к-і(в)^5, если І < 0,

3

1 -£ Аі(і) / Ьі(в)^к-і(в)^в, если і > 0

Дифференцируя обе части последнего равенства по £, приходим к рекуррентному уравнению для функции —::

3 /* о

—й(£)=Е а*(*) / (4)

¿=! ^ — О

Учитывая, что —: является функцией плотности некоторого распределения вероятностей, уравнение (4) следует дополнить условием нормировки

/О

—й(г)Л = 1. (5)

-О

5. Определение предельной плотности —. Будем рассматривать уравнение (4) вместе с условием (5). Указанное уравнение представляет собой линейное интегральное уравнение с вырожденным ядром, решение которого, как известно, может быть сведено к решению некоторой линейной алгебраической системы следующим образом [7, 8]. Для всех *, ] = 1, 2, 3 введем обозначения

/ОО /*О

6; (в)—к (в)^в, 0; = 6*(£)а; (£)Л

-О Л —О

и представим функцию —:(£) в виде

— (*)=£ С^—!) а; (£). (6)

;=!

Умножая обе части равенства (6) на &г(£), а затем интегрируя по всем £, получим

/^(:)

систему рекуррентных уравнений относительно неизвестных С > :

С*0 = Е ^С;(:—!) * =1, 2,3. (7)

;=!

Кроме того, интегрируя обе части (6) по £ и применяя (5), приходим к равенству

Е С(:) = 1. (8)

г=!

Ясно, что имеется взаимно однозначное соответствие между последовательностью функций —:, которая определяется уравнениями (4)—(5), и последовательностью векторов С(:) = (с(:),с2:) ,с3:))т, которая определяется системой (7)—(8), так как каждый шаг, сделанный при переходе от одной системы к другой, обратим.

Исследуем систему (7)—(8). Сначала найдем коэффициенты 0;:

^(11^2 + 9^ + 2^2) ^(2^2 + 7^ + 2^2)

"и = ;—7т- > "12 =

2(^ + ^)(2^ + V)2 ’ (^ + ^)(2^ + ^)(^ + 2^)’

^(11^2 + 10^ + 2^2)

"13 =

2(2^ + ^)2(^ + V)

^(2^2 + 7^ + 2v2)

^(2^2 + 9^ + 1^2)

(^ + V) (^ + 2v)(2^ + V)’ 2(^ + V) (^ + 2v)2

^(2^2 + 10^ + 11v2)

023 =

2(^ + V) (^ + 2v)2

^(11^,2 + 11^ + 2v2) ^(2^,2 + 11^ + 11v2)

31 7777л ! \ о / ! 77 ! "7 7? ^32

2(2^ + v)2(^ + v)(^ + 2v) ’ 2(^ + 2v)2(2^ + v)(^ + V) ’

^(11^,2 + 23^ + 11v2)

0зз — — -

2(^ + 2v)2(2^ + V)2

Покажем, что последовательность векторов, которая определяется системой (7)—(8), сходится. Действительно, с учетом (8) уравнения (7) можно представить в виде

Ср) = £ С;(:—!) + ^, * = 1,2, 3, (9)

;=!

где <г* — некоторая постоянная, 0*; = 0; — а*.

Выберем величину а так, чтобы обеспечить минимум Д2 выражения

А2 = £ £ % -

*=!;=!

После выполнения соответствующих вычислений находим

2 _ К^2^,«6 + 99^ + 27^2 — 83^3 + 27vV + 99v£> + 37v6)

Д22

3(^ + v)2(2^ + V)4 (^ + 2v)4

Как нетрудно проверить, при этом выполняется неравенство Д2 < 1.

Система уравнений (9) задает некоторое отображение трехмерного пространства в себя. Учитывая, что евклидова норма матрицы отображения Д2 < 1, это отображение является сжатием. Следовательно, последовательность векторов С(:) сходится к неподвижной точке отображения при любом начальном векторе С(!).

Неподвижная точка С = (С!, С2, Сз)т является решением системы уравнений

3

с = Е 0; С, * = 1, 2, 3

;=!

при условии

Е С = 1.

г=!

Решая указанную систему, находим С! =

v(2^ + V) (16^6 + 155^ + 539^2 + 850^3 + 652^4 + 240^5 + 32v6)

4(8^8 + 80^ + 321^2 + 690^3 + 880^4 + 690^5 + 321^6 + 80^7 + 8v8) ’

C2 = ц(ц + 2i/)(32/j>6 + 240/rV + 652/x4 v2 + 850/x3z/3 + 539/x2 z/4 + 155 /xz/5 + 16z/6) “ 4 (8/х8 + 80/х7г/ + 321/х6г/2 + 690/х5г/3 + 880/х4И + 690/x3z/5 + 321/x2z/6 + 80/xz/7 + 8г/8) ’ Сз = ^v(^ + 2v)(2^ + v)(8^,4 + 67^3v + 138^2v2 + 67^v3 + 8v4)

4(8^8 + 80^ + 321^2 + 690^3 + 880^4^4 + 690^5 + 321^6 + 80^7 + 8^8)'

В силу взаимно однозначного соответствия между решениями уравнений (4) и (7) и с соответствующими условиями нормировки последовательность функций —^ также сходится к некоторой предельной функции —, которая имеет вид

3

—(г) = Е ^¿(г).

¿=1

Опираясь на представление —^ в виде (6), легко проверить, что последовательность функций —^ сходится к — равномерно.

Наконец ясно, что — является плотностью некоторого распределения вероятностей.

6. Вычисление средней скорости роста. Найдем величину Л по формуле (2). Сначала вычислим интеграл

f G(t — s)^(s)ds — Ci + С2 + С3 +

J —^

2(М — v)

- f /2 -Ci + -G2 + ^ + г/^2 G3^) е“"4 + —С3е-^+^4 - -Crfe-1'*. у /х — v 4 /х(/х + 2г/) J 2ц 2

Тогда функция распределения Ф принимает вид

Ф(*) = F{t) J G(t - s)ijj(s)ds = С\+С2+С3 + (a^Ci - G2 - G3^) e ,2

Cie-Mi —

- (+ 3-c2 + e-‘ - —L—nte-™+

y/x2 — v 4 /x(/x + 2v) J 2(/x — v)

+ (-Ac. + jC2 + •^■;+^+^гсЛ _ ji6.

У^,2 — v2 4 2^(^ + 2v) J 2^

- ^C2te~vt + |G2te~(M+,y)i.

Определим соответствующую функцию плотности:

Ф(*) = м f ——,Ci +с2+ G3^) e-Kf + v( f C\ + —C2 + ^.M + ^2. Gsl e-"4+

^\2(n-v) J 4 /x(/x + 2z/) J

+ — (p + iy) f-J^Cl + -^G2 + ^^G3) e-(^)4+

^ — v y^2 — v2 4(^ + v) 2^(^ + 2v) J

v

Применяя (2), найдем

4/t3 + 8/х2г/ + 11/хг/2 + 5г/3 5/х3 + 11/х2г/+ 8/хг/2 + 4г/3

4 /хг/(/х + г/)2 4 /хг/(/х + г/)2

4/х4 + 14/х3г/ + 21/xV + 14/xz/3 + 4г/4 ^ 2 /хг/(/х + v) (/х + 2г/) (2/х + v) 3

Заметим, что при v — ^ получаем ту же величину Л — 407/(288^), что и в работе [4]. Если v ^ ^, то Л « 5/(4^), что соответствует результату, полученному в [5].

Summary

N. K. Krivulin. Evaluation of the growth rate of the state vector in a generalized linear stochastic system.

A second-order generalized linear stochastic dynamical system is considered. The entries of the system matrix are assumed to be independent and exponentially distributed. In order to evaluate the mean growth rate of the system state vector, a sequence of one-dimensional probability distributions is introduced. Derivation of the limiting density function is reduced to solving an algebraic linear system. The density is used for evaluation of the mean growth rate for the system under study.

Литература

1. Baccelli F., Cohen G., Olsder G. J., Quadrat J.-P. Synchronization and linearity. Chichester: Wiley, 1992.

2. Маслов В. П., Колокольцов В. Н. Идемпотентный анализ и его применение в оптимальном управлении. М.: Физматлит, 1994.

3. Kingman J.F. C. Subadditive ergodic theory // Ann. Probab. 1973. Vol. 1. P. 883-909.

4. Olsder G. J., Resing J. A. C., De Vries R. E., Keane M. S., Hooghiemstra G. Discrete event systems with stochastic processing times // IEEE Trans. Automat. Contr. 1990. Vol. 35. N 3. P. 299-302.

5. Jean-Marie A. Analytical computation of Lyapunov exponents in stochastic event graphs // Performance Evaluation of Parallel and Distributed Systems. Solution Methods: Proc. 3rd QMIPS Workshop. Amsterdam: CWI, 1994. P. 309-341. (CWI Tracts, Vol. 106.)

6. Кривулин Н. К. Скорость роста вектора состояний обобщенной линейной динамической системы со случайной треугольной матрицей // Вестн. С.-Петерб. ун-та. 2005. Сер. 1. Вып. 1. С. 33-38.

7. Михлин С. Г. Лекции по линейным интегральным уравнениям. М.: Физматлит, 1959.

8. Петровский И. Г. Лекции по теории интегральных уравнений. М.: Изд-во Моск. ун-та, 1984.

Статья поступила в редакцию 15 марта 2007 г.