CC BY
43
Intellectual Technologies on Transport. 2015. №2
Метод составления систем однородных дифференциальных уравнений для расчёта вероятностно-временных характеристик нестационарных систем обслуживания
Сергеев С.А.
Петербургский государственный университет путей сообщения императора Александра I
Санкт-Петербург, Россия [email protected]
Аннотация. Предлагается метод построения систем дифференциальных уравнений на основе использования рекурсивного алгоритма генерации матрицы коэффициентов для системы однородных дифференциальных уравнений, описывающих модель нестационарной системы обслуживания. Предлагаемый на его основе алгоритм упрощает построение матриц коэффициентов. С его помощью впервые были реализованы такие сложные модели нестационарных систем обслуживания как сетевая и L-канальная. Также приводится детальный алгоритм построения матрицы коэффициентов для сетевой модели.
Ключевые слова решение нестационарных систем обслуживания, алгоритм генерации матрицы коэффициентов, сетевая модель нестационарной системы осблуживания.
Введение
Отличительной чертой современных аппаратнопрограммных комплексов (АПК) является то, что при создании они, прежде всего, должны быть ориентированы на функционирование не только в нормальных, но ив критических (кризисных) условиях. Для определения возможности реализации всех операций, связанных с технологическим циклом управления, на заданном временном интервале применяют математическое моделирование.
Математической базой является теория массового обслуживания (ТМО). Большинство авторов используют модели ТМО в предположении, что очередь заявок бесконечна, существует стационарный режим, а коэффициент загрузки не превышает единицы [1-4]. Однако наибольший практический и теоретический интерес представляют модели нестационарных систем обслуживания (НСО).
В [5] описан численно аналитический метод, в котором скорость и точность решения обыкновенных дифференциальных уравнений (ОДУ), описывающих НСО, превосходят наиболее распространённый, при решении данного рода задач, метод Рунге-Кутты. В [6] описана его модификация и метод борьбы с погрешностью.
Рассмотрим простейшую одноканальную НСО. На вход НСО последовательно поступает N запросов на об-
работку. Распределения временных интервалов между моментами поступления запросов описываются экспоненциальными законами с интенсивностями {А^А2,AN}, где А; соответствует i-й поступающей заявке. Считаем, что система не имеет потерь. Закон распределения времён обслуживания тоже экспоненциальный с интенсивностями {рх, g2,..., gN}, где g соответствует j-й обслуживаемой заявке.
Состояния системы в каждый момент времени характеризуются числом находящихся в системе запросов i = 0, N и числом уже получивших обслуживание запросов j = 0, N — i. Вероятности пребывания системы в этих состояниях обозначается через P;j(t). Их общее число K = (N + 1)(N + 2)/2. На рис. 1 представлена диаграмма переходов между состояниями системы.
Для определения распределения вероятностей нахождения системы обслуживания в состояниях (i,j) необходимо решить относительно P;j(t) систему ОДУ, каждое из которых выглядит следующим образом:
РуОО = u(i) (Ai+jPi-1J(t) — gj+iPi,j(t)) +
+ u(j)gjPi+ij_i(t) — u(N — i — j)Ai+j+iPiJ(t). (1)
Здесь u(t) — функция Хевисайда, заданная как
и
(0 = {
1,
0,
t > о, t < о.
(2)
Рассмотренная система ОДУ (1) является линейной однородной системой уравнений. Она может быть представлена в матричной форме:
x (t) = Ax(t),
(3)
где x(t) — вектор неизвестных функций размерности K, а A — квадратная матрица.
Интеллектуальные технологии на транспорте. 2015. №2
32
Рис. 1. Диаграмма переходов между состояниями
Когда известно общее уравнение системы ОДУ, описывающей НСО, реализовать алгоритм генерации матрицы А не составляет труда. Но, при добавлении в модель НСО новых характеристик, таких как несколько каналов обслуживания, различный приоритет заявок, очереди к каналам и т.д., общее уравнение системы значительно усложняется. Для примера приведём общее уравнение системы ОДУ, описывающей двухканальную систему обслуживания, приведённую в [7].
НСО имеет N (N=2) каналов обслуживания. На НСО поступает К запросов с интенсивностями {Л1,Лк], интенсивность обслуживания
Р = У] = [
^1,1
Ркл
^1,2 ] Рк,2]
(4)
где s=1,K - номер запроса; j = 1, 2 - номер канала.
Состояния системы в каждый момент времени характеризуются числом обслуженных запросов n (n=0, K), числом находящихся в НСО запросов m (m=0, K-n) и вектором занятости каналов l(l = l1, l2 = [Z7-]), где lj = 0,1. (m, n, 1)-тое уравнение принимает вид:
= Кт) * [ЯДч(Рт-ХпЛ. (О +
8(т- ЮРт-1,п1\ (t) * Лт+п + S(n) *
&Д [Pjl lfc=0+”-2 Рш+Xn-Uj (t) * Pk,j + Pj *
Pm+1,lj (0 * Pn+n-1,j + Wj * Pm+1,n-1,lj (0 * Pm+n,j + wa * Ц=™+п-1 Pm+1,n-1,ij (0 * Pk,j + Si (i1) *
Интеллектуальные технологии на транспорте. 2015. №2
33
Intellectual Technologies on Transport. 2015. №2
$1 O2) *Zk = l+n Pm+1,n-1,lj (t) * VkJ ]} - [s(K -m-
n) * Am+n+1 + S(m) YtjZi S(W * bjj] * Pm,n,i(t), (5)
A
AxT.len, Num X, Num, Num ANum,
Num
X,
(6)
где N=2;
Pj = S1(N + n- lj) * S(N + n-k- 1);
Pn = S±{N + n- lj) * St(N + n-lt- 1); aj = dA(m - 2) * d(li - lj) + dA(m - 1) * dA(lj) * pi; Wj = 5l(li) * @(li) *5(m + n- li);
Wji = Siik) * p(k) * St(m + n- li);
. = (1, если j Ф 1;
1 = [ 2,
если j = 1;
1, если x >0; если x < 0.
1, если x = 0; еслих Ф 0.
т = {0’,
5i (х) = {0
^ = (Ui - ш, 12], если j Ф 1; i llA, l2 - т], если j = 1.
= (ih - (N + n) - к, l2], если j Ф 1; l[l1, l2 - (N + n) - к], если j = 1.
ill = l[li - 1, h], если j Ф 1; j U^, l2 - 1], если j = 1.
pll = ([р + m + n, l2], если j Ф 1; l ([l1, l2 + m + n], если j = 1.
dlll = l[li + к, l2], если j Ф 1; l U^, l2 + к], если j = 1.
Вывести такое уравнение и не допустить в нём ошибки очень трудно. При этом существуют куда более сложные системы, вывод общего уравнения для которых не представляется возможным. Например, сетевая модель, описанная в [5]. Для возможности расчёта вероятностно временных характеристик таких моделей было принято решение разработать метод генерации матрицы А без использования общего уравнения системы ОДУ.
СУТЬ МЕТОДА
Суть метода заключается в том, что необходимо вывести правила перехода из одного состояния в другое и ввести некоторое начальное состояние. Под переходом из одного состояния в другое понимается поступление новой заявки на обслуживание, либо обработка заявки, которая в исходном состоянии обрабатывалась. Начальным состоянием может быть состояние системы, в котором n=0 и m=0. На входе в метод передаём начальное состояние и его номер в списке состояний(№т). Для начального состояния Num=1. Далее метод генерирует одно из возможных состояний по правилам перехода, проверяет было ли сгенерировано данное состояние до этого и если нет, то добавляет его в список состояний (xT). Затем записывает изменения в матрицу А по следующей формуле:
где xT.len - длина списка состояний, сгенерированных на данный момент, х - интенсивность перехода из исходного состояния в новое. Ею может быть интенсивность поступления новой заявки, либо интенсивность обработки очередной заявки. Если полученное состояние уже было сгенерировано и находится в списке состояний под номером Num1, то формула (6) принимает вид:
^Numi, Num X;
ANum, Num ANum, Num
(7)
Если сгенерированное состояние не было получено ранее, то метод начинает рекурсивно генерировать состояния, в которые система может перейти из только что полученного. При этом в рекурсию передаётся Num=xT.len. Если для какого либо состояния сгенерированы все состояния, в которые может перейти система, то этот этап рекурсии завершается. Алгоритм прекращает работу, когда завершатся все этапы рекурсии. Блок-схема алгоритма представлена на рис. 2.
АЛГОРИТМ ДЛЯ ДВУХКАНАЛЬНОЙ СИСТЕМЫ Продемонстрируем детальный алгоритм генерации матрицы А на примере модели двухканальной системы, описанной выше.
Для возможности определения номера заявки, которую канал системы обслуживает в данный момент, в модели НСО вектор I был преобразован в такой что, lj = 0, если j-й канал в данный момент свободен, lj = I, если j-й канал в данный момент обрабатывает i-ю заявку.
Так как состояния системы характеризуются значе-
ниями (n,m,1), то в данной системе переходы между состояниями происходят по следующим правилам:
1) Из состояния Рптр n+m <K, система переходит в состояние Рп1т1-ц, где n1=n,m1=m+1, с интен-
2)
сивностью Ап+т+1. При этом может выполняться одно из условий:
a.
Если 11 = 0, то 111 112 = ^2, 112 = ^2 .
п + 1,112 = 0 п + 2,112 > 0
b. Если 12 = 0 и 11Ф 0, (п + 1,111 = 0
in + 2,111 > 0
то
11 = U
1 =
Из состояния Рп mj, 11 Ф 0, система переходит в
состояние Рп1т1л, n1=n+1,m1=m-1, 1'11 = 0,
112 = 12 с интенсивностью цг1д
Интеллектуальные технологии на транспорте. 2015. №2
34
Intellectual Technologies on Transport. 2015. №2
(Входные параметры)
S - Состояние системы xT - Список всех состояний Num - Номер состояния S
Рекурсивный вызов
3) Из состояния Рпт-[, 12 ^ 0, система может перейти в состояние Pn1m1jz, n1=n+1,m1=m-1, 11г = 1-l, 112 = 0 с интенсивностью |Т;2,2.
В качестве начального состояния было выбрано состояние Рп mj>, n=0, m=0, Z = [0,0].
Блок схема алгоритма, реализующего предлагаемый метод для двухканальной системы с интенсивностью поступления заявок, зависящей от номера заявки, и интен-
сивностью обслуживания, зависящей от номера заявки и номера канала обслуживания, представлена на рисунке 3.
При этом в алгоритме используется дополнительная функция isExist, которая возвращает номер состояния S в списке состояний xT или 0, если состояния в списке нет. Её параметрами являются S - состояние, номер которого ищем, xT - список состояний и Num1 - номер в списке состояний, который требуется найти. Алгоритм данной функции не приводится, т.к. это стандартный алгоритм поиска объекта в списке.
Интеллектуальные технологии на транспорте. 2015. №2
35
Intellectual Technologies on Transport. 2015. №2
Интеллектуальные технологии на транспорте. 2015. №2
36
Intellectual Technologies on Transport. 2015. №2
Рис. 3. Алгоритм для двухканальной модели
Интеллектуальные технологии на транспорте. 2015. №2
37
Intellectual Technologies on Transport. 2015. №2
Как видно из примера выше предложенный алгоритм значительно упрощает процесс составления программы для генерации матрицы коэффициентов А для системы дифференциальных уравнений. Результаты работы модели были проверены другой моделью, где матрица коэффициентов А формировалась с помощью общего уравнения системы ОДУ.
Набольший интерес представляет реализация алгоритма генерации матрицы А для сетевой модели, описанной в [7], так как общее уравнение для неё так и не было выведено.
АЛГОРИТМ ДЛЯ СЕТЕВОЙ МОДЕЛИ
Модель характеризуется следующими параметрами:
1. М - количество типов заявок, причём тип заявки определят относительный приоритет у i-го типа приоритет ниже, чем у (i+1)-ro типа i=(1, М).
2. KTj (i = 1 , М) - количество заявок i-го типа, поступающих в НСО за интервал моделирования.
3. Хц - элемент матрицы интенсивности поступления заявки i-го типа с номером j.
4. L— число канало обслуживания.
5. pу (i = 1,М, j = 1,L) - элемент матрицы
интенсивности обслуживания каналом с номером j заявки i-го типа.
6. фц (i = 1,М, j = 1 ,L) - элемент матрицы
доступности каналов. ф^ = 1 заявке типа I доступен канал с номером j, ф^ = 0 в противоположном случае.
7. n = [щ}(1 = 1, М) - вектор, определяющий
количество заявок i-го типа, находящихся в НСО.
8. т = [mi](i = 1,М) - вектор, определяющий количество заявок каждого типа, уже получивших обслуживание и покинувших НСО.
9. I = [lj}(j = 1,L) - вектор состояний каналов обслуживания, lj = 0, если j-й канал свободен, lj = i, если j-й канал обслуживает заявку типа i.
10. dqt - длина очереди к i-му каналу обслуживания, i = Ц.
11. qt = {qij},j = 1,d({l - вектор очереди к i-му каналу обслуживания, при этом j-я компонента i-го вектора равны номеру типа заявки, стоящей в очереди к i-му каналу обслуживания на j-ом месте.
12. pr = [pri}, i = 1 ,М - i-я компонента вектора определяет приоритет заявок i-го типа.
Состояние НСО в каждый момент времени характеризуются упорядоченным набором векторов:
(п, т, I, dQ, Q) (8)
Алгоритм также использует дополнительную функцию isExist. Кроме того, для возможности гибкого распределения задач между узлами (каналами обслуживания) было введено понятие алгоритма выбора канала обслуживания. Алгоритм выбора канала обслуживания реализован в дополнительной функции choiceChannel. Эта функция имеет 3 параметра: ф - матрица доступности каналов; S -текущее состояние системы; zType - номер типа заявки; channel - переменная, в которую запишется результат работы функции.
Переходы между состояниями системы происходят по следующим правилам:
1. Из состояния Pji^idQQ в котором щ + mt < KTt система может перейти в состояние
р^,т,п,тм с интенсивностью лi,ni+mi+i. при этом номер канала (j), на который попадёт заявка, определяется с помощью функции choiceChannel, если lj > 0, то dq1j = dqj + 1, а i добавляется в qj перед всеми заявками, имеющими меньший приоритет. Если lj = 0, то I1j = i, вектор dq и матрица q остаются без изменения. А n1j = nt + 1. Вектор т остаётся без изменений.
2. Из состояния Pjiffii^ в котором lt >0 система
может перейти в состояние Р^^jldQi qi с интенсивностью щ.д, где n1i = — 1, m1i = ть +
при этом необходимо 'удалить qi}1
Кроме того, в алгоритме используется функция проверки наличия состояния isExist, такая же как и для двухканальной системы, описанной выше. В качестве начального состояния выбрано состояние, в котором 0 обслуженных заявок и 0 поступивших. Блок-схема алгоритма генерации матрицы состояний для сетевой модели приведена на рис. 4.
Интеллектуальные технологии на транспорте. 2015. №2
38
Intellectual Technologies on Transport. 2015. №2
(Входные параметры)
Интеллектуальные технологии на транспорте. 2015. №2
39
Intellectual Technologies on Transport. 2015. №2
Рис. 4. Блок-схема алгоритма генерации матрицы состояний для сетевой модели
Интеллектуальные технологии на транспорте. 2015. №2
40
Intellectual Technologies on Transport. 2015. №2
Данный алгоритм был реализован и интегрирован в программу, описанную в [6]. Кроме того, реализована имитационная сетевая модель. В табл. 1 приведены вероятности поглощающего состояния системы в различные моменты времени, полученные c помощью двух моделей.
Таблица 1
Сравнение результатов моделей
Момент Имитационная Аналитическая
времени модель
i 0.010411 0.0105599357
2 0.224826 0.225337103263
3 0.564251 0.56480951
4 0.776363 0.776569477849
5 0.869929 0.87000275886
СОРТИРОВКА ЭЛЕМЕНТОВ МАТРИЦЫ
Как указано в [5], для системы 3 существует аналитическое решение, если только известны собственные числа матрицы А. Очевидным является тот факт, что в случае треугольного вида матрицы А (для определённости будем считать её нижней треугольной), её собственные числа выписаны в явном виде на диагонали. Таким образом, аналитическое решение системы можно легко найти, если только матрица А - треугольная.
Если пронумеровать состояния по возрастанию числа обработанных запросов, а внутри этих групп по возрастанию числа поступивших запросов то ни одно из состояний полученного списка не будет иметь зависимости от последующих, а матрица А примет треугольный вид. Например, в сетевой модели в одной подгруппе одной группы будет более одного состояния, но нумерация останется верной, так как эти состояния не будут зависеть друг от друга.
На выходе предлагаемый метод выдаёт не треугольную матрицу А и массив xT, в котором хранятся все состояния системы. Нам необходимо отсортировать элементы массива так, как это описано выше. Для этого необходимо воспользоваться любым алгоритмом сортировки. При этом, при перестановке двух состояний с номерами х1 и х2 необходимо поменять местами 2 строки и 2 столбца в матрице А с теми же номерами. Пример такой перестановки приведён в табл. 2 и 3.
Таблица 2
Фрагмент матрицы до перестановки
1 2 х1 х2 х3
1 1-1 2-1 х1-1 х2-1 х3-1
2 1-2 2-2 х1-2 х2-2 х3-2
х1 1-х1 2-х1 х1-х1 х2-х1 х3-х1
х2 1-х2 2-х2 х1-х2 х2-х2 х3-х2
х3 1-х3 2-х3 х1-х3 х2-х3 х3-х3
В табл. 3 приведён фрагмент матрицы состояний А. Первый столбец и первая строка содержат условные номера состояний. В ячейках таблицы содержится информа-
ция вида: x-y, обозначающая интенсивность перехода из состояния х в состояние у.
Таблица 3
Фрагмент матрицы после перестановки
1 2 х2 х1 х3
1 1-1 2-1 х2-1 х1-1 х3-1
2 1-2 2-2 х2-2 х1-2 х3-2
х2 1-х2 2-х2 х2-х2 х1-х2 х3-х2
х1 1-х1 2-х1 х2-х1 х1-х1 х3-х1
х3 1-х3 2-х3 х2-х3 х1-х3 х3-х3
Допустим в процессе сортировки появилась необходимость поменять местами состояния с номерами х1 и х2. Матрица коэффициентов А после такой перестановки представлена в табл. 2. Как видно из табл. 2, зависимости между состояниями не нарушены, что доказывает адекватность предлагаемого метода сортировки.
ЗАКЛЮЧЕНИЕ
Какой бы сложной и многокритериальной не была модель НСО, предложенный подход позволяет синтезировать для неё алгоритм генерации матрицы коэффициентов. Увеличение числа параметров системы ведёт к увеличению числа состояний и, иногда, правил переходов, но сложность вывода этих правил практически не меняется. Таким образом, рекомендуются пользоваться данным методом взамен вывода общего уравнения системы ОДУ.
Литература
1. Osogami T. Analysis of transient queues with semidefinite optimization. / T. Osogami, R. Raymond // Queueing Systems. - 2013. №73. - С. 195-234.
2. Wolff R. W. Little’s law when the average waiting time is infinite. / R. W. Wolff, Y. Yao // Queueing Systems. -2014. №76. - С. 267-281.
3. Sudhesh R., Vijayashree K. V. Stationary and transient analysis of M/M /1 G-queues // Int. J. of Mathematics in Operational Research, 2013. Vol. 5, No 2, pp. 282-299.
4. Sudhesh R., Francis Raj L. Stationary and transient solution of Markovian queues — an alternate approach // Int. J. of Mathematics in Operational Research, 2013. Vol. 5, No 3, pp. 407-421.
5. Бубнов В.П. Алгоритм аналитического расчёта вероятностей состояний нестационарных систем обслуживания. / В.П. Бубнов // Известия Петербургского университета путей сообщения. - 2011. №4. - С. 90-97.
6. Бубнов В.П. Особенности программной реализации численно-аналитического метода расчёта моделей нестационарных систем обслуживания. / В.П. Бубнов,
A. С. Еремин, С.А. Сергеев // Труды СПИИРАН. - 2015. №1. - С. 218-232.
7. Бубнов В.П. Разработка динамических моделей нестационарных систем обслуживания. / В.П. Бубнов,
B. И. Сафонов. - г. Санкт-Петербург, 1999, 65 с.
Интеллектуальные технологии на транспорте. 2015. №2
41
Intellectual Technologies on Transport. 2015. №2
The Method of Construction of Systems of Homogeneous Differential Equations for Calculating the Probability-Time Characteristics of non-Stationary Service Systems
Sergeev S.A.
Petersburg State Transport University Saint-Petersburg, Russia [email protected]
Abstract. It is proposed the method of construction of systems of differential equations based on the use of a recursive algorithm for generating the coefficient matrix for the homogeneous system of differential equations describing the model of non-stationary service systems. The proposed approach simplifies the construction of the matrix coefficients. With it was first implemented such complex models of non-stationary systems as network and L-channel models. The article also provides a detailed algorithm for the network model.
Keywords: solution of non-stationary systems, service, algorithm for generating the coefficient matrix, network model non-stationary system inspection completed.
REFERENCES
1. Osogami T. Analysis of transient queues with semidefinite optimization. / T. Osogami, R. Raymond // Queueing Systems. 2013. № 73. pp. 195-234.
2. Wolff R. W. Little’s law when the average waiting time is infinite. / R. W. Wolff, Y. Yao // Queueing Systems. -2014. № 76. pp. 267-281.
3. Sudhesh R., Vijayashree K. V. Stationary and transient analysis of M/M/1 G-queues
// Int. J. of Mathematics in Operational Research, 2013. Vol. 5, No 2, pp. 282-299.
4. Sudhesh R., Francis Raj L. Stationary and transient solution of Markovian queues — an alternate approach // Int. J. of Mathematics in Operational Research, 2013. Vol. 5, No 3, pp. 407-421.
5. Bubnov V.P. Algoritm analiticheskogo raschyota veroyatnostej sostoyanij nestacionarnyh system obsluzhivani-ya: Izvestiya Peterburgskogo universiteta putej soobshcheni-ya. [Algorithm of analytical calculation of non-stationary state probabilities service systems: News from the St. Petersburg University of communication.] / V.P. Bubnov. 2011. № 4. 9097 p.
6. Bubnov V.P. Osobennosti programmnoj realizacii chislenno analiticheskogo metoda raschyota modelej nes-tacionranyh sistem obsluzhivaniya: Trudy SPIIRAN. [Features of the software implementation of numerical-analytical method of calculation models non-stationary service systems: SPIIRAS Proceedings.] / V.P. Bubnov, A.S. Eremin, S.A. Sergeev. 2015. №1. pp. 218-232.
7. Bubnov V.P. Razrabotka dinamicheskih modelej nestacionarnyh system obsluzhivaniya. [Developing dynamic modeling of non-stationary systems.] / V.P. Bubnov, V.I. Safonov. - Saint-Petersburg, 1999, 65 p.
Интеллектуальные технологии на транспорте. 2015. №2
42
