Имитационное моделирование систем и сетей массового обслуживания средствами приложения к пакету Mathcad

УДК 004.94

Е. Б. Захарикова, П. П. Макарычев

ИМИТАЦИОННОЕ МОДЕЛИРОВАНИЕ СИСТЕМ И СЕТЕЙ МАССОВОГО ОБСЛУЖИВАНИЯ СРЕДСТВАМИ ПРИЛОЖЕНИЯ К ПАКЕТУ MATHCAD

Аннотация. Рассматриваются особенности разработанного пакета прикладных программ имитационного моделирования, реализованного на языке С++ и встроенного в Mathcad. Проводится анализ адекватности, чувствительности и устойчивости конкретной модели.

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

Abstract. The article considers the features of the developed software package of simulation, implemented in C++ and built into Mathcad. The authors analyze adequacy, sensitivity and stability of the model.

Key words: queuing system, queuing network, analytical model, simulation model, Markov process, labeled transition graph of a random process, system of differential equations of Kolmogorov, Wilcoxon test.

Введение

В процессе создания и исследования сложных систем важное место занимает имитационное моделирование. Реализация имитационного моделирования может быть осуществлена с использованием специализированных проблемно-ориентированных программных средств (Arena, AnyLogic, GPSS, Vissim, ExtendSim, AutoMod, Promodel и др.) и универсальных математических пакетов (Mathcad, Matlab, Maple, Mathematica и др.).

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

1. Описание структуры и принципа работы пакета прикладных программ имитационного моделирования

Модель структуры СеМО задается в виде ориентированного графа, вершины которого представляют множество возможных узлов обслуживания. В качестве основных компонентов сети рассматриваются:

- генератор, предназначенный для моделирования входного потока заявок по заданным законам распределения (экспоненциальный, нормальный, детерминированный или эмпирический). Количество генераторов заявок не ограничено;

- узел обслуживания, осуществляющий обслуживание заявок по заданным законам распределения (экспоненциальный, нормальный или детерминированный). Количество узлов обслуживания не ограничено;

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

Расчет статистических характеристик выполняется на основе накопления данных и многократного повторения эксперимента.

Данный комплекс программ имеет модульную структуру и встроен в Mathcad с помощью подключения внешней библиотеки динамической компоновки DLL (Dynamic Link Library), написанной на языке C++. При этом пользователь автоматизированной системы моделирования располагает возможностями абсолютного контроля над своей моделью, может варьировать по желанию любой параметр и судить о поведении модели по наблюдаемым результатам.

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

2. Анализ адекватности имитационной модели

Оценка адекватности модели была выполнена методом предельных точек путем сравнения результатов, полученных при аналитическом и имитационном моделировании [1]. В качестве примера была рассмотрена двухузловая СМО без буфера с неоднородным потоком заявок от двух источников с интенсивностями Х1 = 0,4 и = 0,3 соответственно. Интенсивности обслуживания ц1 = 0,5 и ц2 = 0,6 зависят от типа заявки и не зависят от узла. Концептуальная модель данной системы в международном стандарте моделирования IDEF3 представлена на рис. 1.

Рис. 1. Концептуальная модель

Сформулируем предположения и допущения для аналитического моделирования:

- поступающие в систему заявки двух типов образуют простейшие потоки;

- длительности обслуживания заявок каждого типа распределены по экспоненциальному закону;

- дисциплина обслуживания FIFO: заявка, поступившая в систему, принимается на обслуживание, если есть хотя бы один свободный узел. Если заявка застала свободными несколько узлов, то она направляется в один из них случайным образом.

Случайный процесс, протекающий в системе с дискретными состояниями и непрерывным временем, является марковским.

Для описания состояний марковского процесса используется распределение заявок между узлами. Закодируем состояния следующим образом: (У1,У2), где У1, У2 = {0, 1, 2} - состояние обслуживающих приборов, задаваемое типом заявки, находящейся на обслуживании («0» - узел свободен; «1» или «2» - на обслуживании в приборе находится заявка типа 1 или 2 соответственно).

При выбранном способе кодирования система может находиться в следующих состояниях:

Е0:(0,0) - в системе нет ни одной заявки;

Ej:(1,0) - на обслуживании в узле 1 находится заявка типа 1;

E2:(2,0) - на обслуживании в узле 1 находится заявка типа 2;

E3:(1,1) - на обслуживании в узле 1 находится заявка типа 1 и на обслуживании в узле 2 находится заявка типа 1;

E4:(1,2) - на обслуживании в узле 1 находится заявка типа 1 и на обслуживании в узле 2 находится заявка типа 2;

E5:(0,2) - на обслуживании в узле 2 находится заявка типа 2;

E6:(0,1) - на обслуживании в узле 2 находится заявка типа 1;

E7:(2,2) - на обслуживании в узле 1 находится заявка типа 2 и на обслуживании в узле 2 находится заявка типа 2;

E8:(2,1) - на обслуживании в узле 1 находится заявка типа 2 и на обслуживании в приборе 2 находится заявка типа 1.

Размеченный граф переходов случайного процесса представлен на рис. 2.

Рис. 2. Размеченный граф переходов случайного процесса

В каждый момент времени может произойти только одно событие (или поступление заявки какого-либо типа, или завершение обслуживания заявки, находящейся в узле), поскольку вероятность появления двух и более событий в один и тот же момент времени равна нулю.

Для вероятностей состояний р0(0, рКО, •••, Р8(0 можно по графу переходов составить систему линейных дифференциальных уравнений Колмогорова [2]:

Р°( ) =М-1 р1({) + М-2р2(() + М-1 р6(() + М-2р5(() -(^1 +^2)р0((),

Ш

р^( ) = ^1р0(^) + М-2р4(^) + М-1 р3()-(^1 + ^2 +М-1)р1(^),

Ш

) = ^2р0(^) + М-1р8(^) + М2р7(()-(^1 + ^2 + М2)р2((X

Ш

) = ^1 р1 (() + ^1 рб (() - 2М-1 рз ((),

Ш

" ) =^2р1

Ш

Р5' ) =М1 Р4(() + М2Р7(()-(^1 + ^2 +^2)Р5(^X

Ш

) =М1 р3() + М2р8() -(^1 + ^2 +М1)рб(),

Ш

) = ^2р5(() + ^2р2(() -2М-2р7((X

Ш

Р8' ) =^2рб(() + ^1 р2(()-(М-1 + М-2)р8(^).

Ш

Начальные условия для интегрирования отражают состояние системы в начальный момент: ро(0 = 1; рг(0 = 0, г = 1, ..., 8.

Было произведено решение данной системы в Mathcad с получением наглядной графической интерпретации, представленное на рис. 3.

Шрг ^)

При предельном режиме ^ ------= 0; г = 0.. .8, полученная систе-

Ш

ма принимает вид системы линейных алгебраических уравнений. Получены следующие значения предельных вероятностей: р0 = 0,318, р! = 0,181, р2 =

0,116, р3 = 0,102, р4 = 0,0б5, р5 = 0,043, рб = 0,073, р7 = 0,04, р8 = 0,062.

Предельная вероятность того, что первый узел занят, равна

р1 = Р1 + Р2 + рз + Р4 + Р7 + Р8, р1 = 0,565.

Предельная вероятность того, что второй узел занят, равна

pH = рз + Р4 + Р5 + Рб + Р7 + р8, р2г = 0,385.

Предельная вероятность потерь равна

Р

С помощью программных средств имитационной модели также были рассчитаны значения вероятностей занятости первого и второго узлов p1p, p2p и вероятность потерь pp [3]. При общем времени моделирования T = 1000 были вычислены времена работы каждого узла с повторением вычислений n = 10 раз при t = 0, 10, ..., T. Листинг главной программы, управляющей процессом моделирования, приведен на рис. 4, где Ns - количество источников; Nk - количество узлов обслуживания; distInStr, DistCan - законы распределения времен генерирования и обслуживания; distPar, ParCan - параметры распределения. Подпрограмма InputStream генерирует входной поток заявок от каждого источника, подпрограмма CombineNStream формирует единый поток заявок, AutoServ инициирует обслуживание.

Графики теоретической и практической вероятностей потерь pt(t), pp(t) и занятости узлов p1t(t), p1p(t), p2t(t), p2p(t) представлены на рис. 5, 6.

Из графиков, представленных на рис. 5 и рис. 6, следует, что в установившемся режиме относительные погрешности вероятности потерь и занятости узлов не превышают 3 и 1,5 % соответственно.

Сравнение результатов моделирования системы приведено в табл. 1.

Таким образом, сопоставление полученных теоретических и практических результатов исследования показало адекватность рассматриваемого пакета прикладных программ имитационного моделирования.

Доступность и наглядность работы системы достигаются визуализацией.

На рис. 7 представлены графики, отражающие функционирование системы, где f(t) и h(t) - распределение во времени заявок первого и второго типов; w(t) и v(t) - распределение во времени занятости первого и второго узлов; w(t) + v(t) - распределение во времени суммарной занятости узлов; p(t) - распределение потерь.

Т := 1000 Ns := 2 Nk := 2 n := 10

distlnStr :=

DistJCan :=

"exp"

/exp"

"exp"

/exp"

distPar :=

PaiCan :=

0.4 0

,0.3 0,

'0.5 O' L0.6 0

result :=

for j € 0..11- 1

x InputStream

s sum

comstr s if Ns = 1

CombineN3tream(s) otherwise

service AutoServ

result. <— service l

result

Рис. 4. Листинг главной программы

Рис. 5. Графики вероятностей потерь

а)

а)

Рис. 6. Графики вероятностей занятости узлов

t 6)

б)

Рис. 6. Окончание

Таблица 1

Сравнение результатов моделирования

Т = 1000 11 = 0,4 Х2 = 0,3 Ці = 0,5 Д2 = 0,6 Параметр Аналитическое моделирование Имитационное моделирование Относительная погрешность, %

Вероятность занятости первого канала 0,565 0,567 0,35

Вероятность занятости второго канала 0,385 0,39 1,3

Среднее время обслуживания заявки первого типа 2 1,986 0,7

Среднее время обслуживания заявки второго типа 1,667 1,633 2

Вероятность потерь 0,269 0,276 2,6

Наряду со стандартной дисциплиной обслуживания многоузловых СМО разработана экономичная дисциплина обслуживания, в которой при наличии нескольких свободных узлов заявка направляется в узел, который к этому времени простаивал больше всего. При этом достигается примерно одинаковая загрузка всех узлов и уменьшается вероятность потерь. Так, для представленного выше примера вероятности занятости первого и второго узлов равны pip = 0,468, p2p = 0,463, вероятность потерьpp = 0,235.

3. Анализ чувствительности имитационной модели

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

ния выходных характеристик математическом и имитационном моделей при незначительном изменении входных характеристик (табл. 2).

W0О+уГО

10 30 30 4 Я 60 ТО ВО PQ 100

1 »

в)

Рис. 7. Графики распределения заявок и потерь в компонентах СМО

Таким образом, относительная погрешность выходных параметров имитационной системы не превышает 3,5 %, что говорит об удовлетворительной чувствительности рассматриваемого пакета прикладных программ имитационного моделирования.

4. Анализ устойчивости имитационной модели

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

какого-либо подмножества генеральной совокупности (т.е. выборки). Устойчивость результатов моделирования можно рассматривать как признак, подлежащий оценке. Для проверки гипотезы об устойчивости результатов может быть использован критерий Уилкоксона.

Таблица 2

Результаты анализа модели на чувствительность

Входные данные Вероятность занятости Среднее время обслуживания заявки Вероятность потерь

первого канал второго канала первого типа второго типа

1 Математическая модель А1 = 0,4 А2 = 0,3 Ц1 = 0,6 Д2 = 0,7 0,523 0,329 1,667 1,429 0,223

Имитационная модель 0,527 0,324 1,679 1,431 0,22

Относительная погрешность 0,7 1,5 0,7 0,1 1,3

2 Математическая модель А1 = 0,4 А2 = 0,3 № = 1 Д2 = 1,2 0,394 0,182 1 0,833 0,113

Имитационная модель 0,398 0,186 0,99 0,85 0,117

Относительная погрешность 1 2,2 1 2 3,5

3 Математическая модель А1 = 0,3 А2 = 0,2 Ц1 = 0,5 Д2 = 0,6 0,483 0,279 2 1,667 0,184

Имитационная модель 0,494 0,283 2,061 1,672 0,188

Относительная погрешность 2,3 1,4 3 0,3 2,2

Критерий Уилкоксона служит для проверки того, относятся ли две выборки к одной и той же генеральной совокупности.

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

Проверка указанной гипотезы проводится в следующем порядке:

1. Формируются две выборки: вероятности занятости первого канала pz11 при АЛ = 0,4 и вероятности занятости первого канала pz12 при АЛ = 0,1 при сохранении значений остальных входных характеристик.

2. Относительно законов распределения pz11 и pz 12 никаких предположений не делается.

3. Значения обеих выборок упорядочиваются по возрастанию (pz11s, pz12s).

4. В случае pz12si < pz11si, i = 0.. .100, пара значений образует инверсию.

5. Подсчитывается полное число инверсий V; V = 100.

6. Если гипотеза верна, то V не должно сильно отклоняться от своего математического ожидания М = т, где т - размер выборки, т = 101.

Очевидно, что в данном случае согласно критерию Уилкоксона выборки принадлежат одной и той же генеральной совокупности, что подтверждает достаточную устойчивость системы. Принадлежность выборок к одной и той же генеральной совокупности иллюстрируют графики, представленные на рис. 8.

Рис. 8. Визуализация анализа модели на устойчивость Заключение

Таким образом, исследование полученных теоретических и практических результатов СМО показало удовлетворительную работу данного пакета прикладных программ имитационного моделирования. Далее планируется проведение исследований с более сложными СеМО.

Список литературы

1. Шеннон, Р. Имитационное моделирование систем. Искусство и наука / Р. Шеннон. - М. : Мир, 1978. - 418 с.

2. Клейнрок, Л. Теория массового обслуживания / Л. Клейнрок. - М. : Машиностроение, 1979. - 432 с.

3. Кельтон, В. Имитационное моделирование / В. Кельтон, А. Лоу. - СПб. : Питер, 2004. - 847 с.

4. Алгазинов, Э. К. Анализ и компьютерное моделирование информационных процессов и систем / Э. К. Алгазинов, А. А. Сирота. - М. : Диалог-МИФИ, 2009. -416 с.

Захарикова Елена Борисовна аспирант, Пензенский государственный университет

E-mail: makpp@yandex.ru

Макарычев Петр Петрович

доктор технических наук, профессор, заведующий кафедрой математического обеспечения и применения ЭВМ, Пензенский государственный университет

E-mail: makpp@yandex.ru

Zakharikova Elena Borisovna Postgraduate student,

Penza State University

Makarychev Petr Petrovich Doctor of engineering sciences, professor, head of sub-department of computer application and software,

Penza State University

УДК 004.94 Захарикова, Е. Б.

Имитационное моделирование систем и сетей массового обслуживания средствами приложения к пакету Mathcad / Е. Б. Захарикова, П. П. Макарычев // Известия высших учебных заведений. Поволжский регион. Технические науки. - 2012. - № 3 (23). - С. 25-35.