Численный анализ системы массового обслуживания с гиперэкспоненциальным «Охлаждением»

ВЕСТНИК ТОМСКОГО ГОСУДАРСТВЕННОГО УНИВЕРСИТЕТА

2016 Управление, вычислительная техника и информатика № 4

УДК 519.872

DOI: 10.17223/19988605/37/4

В.А. Лохвицкий, А.В. Уланов

ЧИСЛЕННЫЙ АНАЛИЗ СИСТЕМЫ МАССОВОГО ОБСЛУЖИВАНИЯ С ГИПЕРЭКСПОНЕНЦИАЛЬНЫМ «ОХЛАЖДЕНИЕМ»

Рассматривается многоканальная система массового обслуживания, в которой после ухода последней заявки периода непрерывной занятости на некоторое время наступает режим «охлаждения». Приходящие за время охлаждения заявки накапливаются и ожидают, пока система не перейдет в режим обслуживания. Получено распределение числа заявок и преобразование Лапласа-Стилтьеса времени ожидания для системы с простейшим входящим потоком, экспоненциальным обслуживанием и гиперэкспоненциальным (с произвольным типом параметров) охлаждением.

Ключевые слова: системы массового обслуживания с «охлаждением»; немарковские системы массового обслуживания; гиперэкспоненциальное распределение; численные методы.

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

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

Возможен и другой случай, когда после полного освобождения системы требуется проведение мероприятий, связанных с восстановлением работоспособности, проведением технического обслуживания, перерывами в работе и т.д. В этом случае уместно говорить, что СМО переходит в режим «охлаждения». Начало обслуживания вновь прибывшей заявки не начнется, пока не завершатся все операции охлаждения. В отличие от разогрева, процесс охлаждения системы не зависит от прибытия первой заявки периода занятости. Если система успеет «охладиться» до прихода заявки, обслуживание начнется без дополнительных задержек.

Из работ по исследованию систем с охлаждением авторам встречалась лишь статья [3], в которой рассматривалась полностью марковская система M/M/M/n. Ниже представлены модель и методика расчета вероятностно-временных характеристик системы с простейшим входящим потоком, марковским обслуживанием и немарковским охлаждением. Время охлаждения аппроксимируется гиперэкспоненциальным распределением второго порядка (Н2) с возможностью комплексного типа параметров, что обеспечивает достаточную для инженерных расчетов точность [4-8]. Модель системы M/M/H2/n представляет собой дискретный марковский процесс с непрерывным временем.

1. Модель системы MIMIH2ln

Рассмотрим модель многоканальной СМО с «охлаждением» типа M/M/H2/n. Напомним, что гиперэкспоненциальное распределение второго порядка относится к распределениям фазового типа и предполагает выбор случайным процессом одной из двух альтернативных фаз [2]. С вероятностью y про-

цесс охлаждения попадает в первую фазу и задерживается в ней случайное время, распределенное по экспоненциальному закону с параметром дС . С вероятностью у2 = 1 - у процесс попадает во вторую фазу, где экспоненциальная задержка имеет параметр д2. Диаграмма переходов между состояниями марковского процесса, описывающего систему Ы/Ы/И2/п, представлена на рис. 1.

&

К

«

о Л

н о К о т и о

т «

СО

о ч

о К

п-1

х\

X

X

(п-1)д\

'X

д 1

п+1

-я

пд{

, д Iе

пд(

\ д 1

-(

д 2

X

X

X

X

к с д 2 / IX

к К

д 2 С

Рис 1. Диаграмма переходов для системы с охлаждением Ы1М1И21п

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

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

2. Методика расчета вероятно-временных характеристик

Обозначим через 8/ множество всех возможных микросостояний системы, при которых на обслуживании находится ровно / заявок, а через о/ - количество элементов в </. Далее, в соответствии с диаграммой переходов построим матрицы интенсивностей инфинитезимальных переходов:

- АДо/ХО/н] - в <8+1 (по прибытию заявок);

- В/ [о/Хо/л] - в / (по завершению обслуживания);

- СДо/Хо/] - в 8/ (по завершению охлаждения);

- £>/[о/Хо/] - ухода из микросостояний/-го яруса (диагональная матрица).

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

Для СМО М/М/И2/п матрицы переходов будут иметь следующий вид:

"X 0 0"

Aj =

0 X 0 0 0 X

j = 0, N -1;

" 0 0 0 " "0 0 0" "0 0 0"

B = ЦУ1 0 ЦУ2 , B2 = 0 2ц, 0 ... B = 0 пц 0

0 0 0 0 0 0 0 0 0

, j = n, N;

=

M1

0

c

Ц 2

D =

X + 0 0

0

, j = 0, N.

0

X + шт(п,/)ц 0

о х +

Очевидно, матрицы А и С имеют одинаковый вид для всех ярусов диаграммы, а матрицы В и В зависят от номера яруса, но при/ > п стабилизируются.

Введем векторы-строки у;- = [у;д, у;;2, у;;3] нахождения СМО в микросостояниях /-го яруса. Теперь запишем векторно-матричные уравнения баланса переходов между состояниями [2], указанными на диаграмме:

У 0 В0 = У0С0 + У1В1,

у/В/ = У/,1 А/1 + У/С/ + У/+1В/+1, /' =1, я,

где N - число обсчитываемых ярусов диаграммы.

Для решения данной системы уравнений воспользуемся итерационным методом, предложенным в [8] и доработанным в [9]. Направление прогонки при расчёте определим сверху вниз - в сторону увеличения номера яруса. В результате расчета получим отношения смежных вероятностей числа заявок в

системе x,

- Pj + i/^ j = 1, N.

После окончания итераций, используя значения X/, осуществим переход к вероятностям состояний системы по следующему алгоритму:

- вероятность свободного состояния р0 положить равной 1;

- рассчитать р/ + 1 = р/ х/, / = 1, N;

- вычислить сумму 8 = 2 Я=0 Р /;

нормировать полученные значения: Pj = Pj /S, j = 0, N.

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

Время ожидания вновь прибывшей заявки определяется состоянием системы непосредственно перед ее прибытием. В соответствии с теоремой PASTA (Poisson Arrival See Time Average) [10, 11] рас-

пределение числа заявок перед прибытием очередной совпадает со стационарным. Перед прибытием очередной заявки система находится в одном из трех микросостояний - охлаждение 1-го или 2-го типов или обслуживание.

Пусть перед прибытием очередной заявки в системе уже находилось у заявок. Вероятности соответствующих микросостояний записаны в векторе-строке у у = [ууд, у;;2, у;;3]. Если система находилась в первом микросостоянии у-го яруса (см. рис. 1), вновь прибывшая заявка будет ожидать завершения охлаждения первого типа плюс ( - п + 1) продвижений очереди. Преобразование Лапласа-Стилтьеса (ПЛС) соответствующего времени

а у,1( =

С пц ^

у-п+1

ц +5

^пЦ + Я^

В третьем микросостоянии происходит охлаждение второго типа. ПЛС соответствующего ожида-

ния

а 7,З(5) = -

Ц2

(

пц

\

у-п+1

пц + 5

ц2 + 5

Второе микросостояние соответствует режиму обслуживания. ПЛС времени ожидания вновь прибывшей заявки составит

Ю у ,2(5) =

С \

пд

ч пд + 5 у

у-п+1

Введем вектор-столбец ю/я) = [юу-д^), ю;;2(5), ю;;з(5)]т. Тогда ПЛС искомого времени ожидания

начала обслуживания

N-1

Ю ( я) = I У у Ю у ( я).

}=п

Для получения начальных моментов времени ожидания можно выполнить дифференцирование в нуле выражения (1) либо заменить все входящие в него ПЛС наборами начальных моментов, а их произведения - свертками в моментах [2].

3. Результаты расчетов

Представленные диаграмму переходов и уравнения баланса переходов между ее микросостояниями принято считать справедливыми для системы Ы/Ы/И2/п лишь в том случае, когда параметры гиперэкспоненты вещественны. Это соответствует коэффициенту вариации времени охлаждения и > 1. Когда же коэффициент вариации времени охлаждения и < 1, подбираемые по методу моментов параметры ^-распределения уходят либо в комплексную, либо в «парадоксальную» область [8]. Строго говоря, в этом случае мы имеем дело уже не с гиперэкспонентой, а с некоторой трехпараметрической функцией

Я (х) = Ле-д'х + у2е-д2х; у + У2 = 1, параметры которой могут принимать комплексные значения. Очевидно, ^-распределение является частным случаем функции Я(х) с вещественными параметрами.

Соответственно, при использовании Я-аппроксимации с комплексными параметрами элементы матриц интенсивностей и параметры уравнений баланса также принимают комплексные значения. Тем не менее многочисленные вычислительные эксперименты [2, 6-8] показывают, что при последующем решении уравнений баланса указанная патология проявляется лишь в промежуточных результатах - вероятностях «фиктивных» микросостояний, на которые расщепляются «фактические» состояния СМО. На этапе суммирования вероятностей микросостояний комплексные части аннигилируются и результат расчета - распределение числа заявок в системе - становится вещественным.

Выполним расчет системы Ы/Ы/И2/п на имитационной модели и численным методом через #2-аппроксимацию для следующего набора исходных данных: число каналов обслуживания п = 3, средняя интенсивность входящего потока X = 2,5, среднее время обслуживания заявок Ъ = 1, среднее время

охлаждения системы С1 = 2. Параметры ^-распределения подбирались по трем моментам различных исходных распределений длительности охлаждения:

- вырожденного В (коэффициент вариации и = 0);

- равномерного и на интервале [0; 4] (и ~ 0,577);

- экспоненциального Ы (и = 1);

- гамма Г с параметром формы 0,5 (и ~ 1,41);

- Вейбулла Ж с параметром формы 0,46 (и = 2,5).

Результаты расчета распределения числа заявок в системе для В, Ы и Ж распределений приведены на рис. 2. Штриховкой показаны графики, полученные с помощью имитационной модели, сплошными линиями - на основе численного метода расчета.

]

Рис 2. Распределение числа заявок в системе с охлаждением Ы/М/И2/п

Из графиков видно согласие результатов даже в области комплексных и парадоксальных параметров гиперэкспоненты. Расстояние Колмогорова между результатами, полученными методом имитационного моделирования и с помощью Я-аппроксимации В, и, Ы, Г и Ж распределений времени охлаждения, составило {0,006; 0,005; 0,006; 0,004; 0,043} соответственно, что говорит о приемлемой точности Я-аппроксимации. Заметим, что для В, и распределений параметры Я-функции принимали комплексные значения, что, тем не менее, не повлияло на корректность конечного результата.

Увеличение расстояния Колмогорова для Ж-распределения с параметром формы 0,46 может быть обусловлено тем, что оно относится к классу распределений с «толстым хвостом» и, следовательно, описываемая им случайная величина имеет очень большой разброс. В этой связи возрастают требования к качеству используемых в имитационной модели программных датчиков случайных чисел.

В таблице представлены результаты расчета начальных моментов распределения времени ожидания {^}, / = 1, 3 , полученные численно (числ.) и с помощью имитационной модели (ИМ).

Начальные моменты времени ожидания

Е(г) м>2

Числ. ИМ Числ. ИМ Числ. ИМ

В 1,788 1,709 6,674 6,495 37,19 34,16

и 1,971 1 ,870 7,864 7,230 40,46 37,22

Ы 2,272 2,225 1,035 1,035 62,88 64,02

Г0,5 2,582 2,494 1,406 1,356 108,0 103,9

Ж 2,600 2,802 2,586 2,476 767,0 687,7

Из представленных результатов следует, что с увеличением коэффициента вариации распределения времени «охлаждения» среднее время ожидания возрастает. Сходство результатов, полученных численно и с помощью имитационного моделирования, позволяет сделать вывод о корректности предложенной методики численного расчета вероятностно-временных характеристик СМО с «охлаждением».

Заключение

Представленные численные результаты свидетельствуют о существенном влиянии коэффициента вариации времени охлаждения на вероятностно-временные характеристики СМО с «охлаждением». Используемая для аппроксимации немарковской компоненты СМО ^-функция позволяет с приемлемой точностью проводить анализ систем массового обслуживания с произвольным коэффициентом вариации распределения времени охлаждения.

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

- восстановлением работоспособности системы;

- проведением технического обслуживания;

- восстановлением ресурсов, расходуемых при обслуживании;

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

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

ЛИТЕРАТУРА

1. Гиндин С.И., Хомоненко А.Д., Ададуров С.Е. Численный расчет многоканальной системы массового обслужи-

вания с рекуррентным входящим потоком и «разогревом» // Известия Петербургского университета путей сообщения. 2013. № 4 (37). С. 92-101.

2. Рыжиков Ю.И. Алгоритмический подход к задачам массового обслуживания. СПб. : ВКА им. А.Ф. Можайско-

го, 2013. 496 с.

3. Хомоненко А. Д., Халиль М. М., Гиндин С. И. Моделирование облачных вычислений с использованием многока-

нальной системы массового обслуживания с «охлаждением» // Материалы XIX Международной конференции по мягким вычислениям и измерениям (SCM-2016). 25-27 мая 2016 г. СПб. : Изд-во СПбГЭТУ «ЛЭТИ», 2016. Т. 1. С. 247-251.

4. Кендалл М. Дж., Стьюарт А. Теория распределений: пер. с англ. М. : Наука, 1966. 587 с.

5. Cox D.R. A Use of Complex Probabilities in the Theory of Stochastic Process // Proc. of the Cambridge Phil. Soc. 1955.

P. 313-323.

6. Рыжиков Ю.И., Уланов А.В. Расчет гиперэкспоненциальной системы M/H2/n-H2 с заявками, нетерпеливыми в

очереди // Вестник Томского государственного университета. Управление, вычислительная техника и информатика. 2014. № 2 (27). C. 47-53.

7. Рыжиков Ю.И., Лохвицкий В. А. Расчёт узла сети с отрицательными заявками // Сборник трудов по материалам

российской конференции «Информационные технологии в управлении» (ИТУ-2012). СПб., 2012. С. 354-359.

8. Рыжиков Ю.И., Уланов А.В. Применение гиперэкспоненциальной аппроксимации в задачах расчета немарков-

ских систем массового обслуживания // Вестник Томского государственного университета. Управление, вычислительная техника и информатика. 2016. № 3 (36). С. 60-65.

9. Takahashi Y., Takami Y. A Numerical Method for the Steady-State Probabilities of a GI/G/c Queuing System in a Gen-

eral Class // J. of the Operat. Res. Soc. of Japan. 1976. V. 19, No. 2. Р. 147-157.

10. Рыжиков Ю. И., Хомоненко А. Д. Итерационный метод расчёта многоканальных систем с произвольным законом обслуживания // Проблемы управления и теории информации. 1980. Т. 9, № 3. С. 32-38.

11. Serfozo R. Basics of Applied Stochastic Processes. Berlin : Springer, 2009. 443 p.

Лохвицкий Владимир Александрович, канд. техн. наук. E-mail: vovan296@mail.ru Уланов Александр Викторович, канд. техн. наук. E-mail: ulanov246@rambler.ru Военно-космическая академия имени А.Ф. Можайского

Поступила в редакцию 23 июня 2016 г.

Lokhvitsky Vladimir A., Ulanov Alexander V. (Mozhaisky Military Space Academy, Russian Federation, St. Petersburg).

The numerical analyses of queuing system with hyperexponential distribution of cooling time.

Keywords: queuing systems with cooling; non-markovian queuing systems; hyperexponential distribution; numerical

methods.

DOI: 10.17223/19988605/37/4

The numerical analysis of QS with cooling is considered. Cooling is understood as the temporary delays coming after the period of continuous employment. In practice such delays arise when after service of the last request in system carrying out the actions connected with maintenance, carrying out maintenance, breaks in work, etc. is required. In this case pertinently to say that QS passes into the cooling mode. The beginning of service of again arrived request won't begin, all operations of cooling won't come to the end yet. Unlike a warming up, process of cooling of system doesn't depend on arrival of the first request of the period of employment. If the system manages to be cooled before arrival of this request, service will begin without additional delays.

In this paper is presented the technique of the numerical analysis of multichannel QS with the Poisson arrival, exponential service and not-exponential distribution of cooling time.

For approximation of cooling time it is offered to use the function R(x) which is formally coinciding with function of the hyperex-ponential distribution of the second order (H2), but allowing unlike hyper exhibitors complex type of parameters

R(x) = y^ + y2e-^x; y + y2 = 1. Considering that prototype of function R(x) - the hyperexponential H2 distribution belongs to distributions of phase type, conditions of non-Markov QS and transitions between them can be presented in the form of discrete Markov process with continuous time. Possibilities of R-approximation in that case, when coefficient of a variation of initial distribution u<1 are shown in article. Function parameters thus accept complex values. Nevertheless, at calculation of QS with application of R-approximation in the field of complex values this pathology is shown only in intermediate results - probabilities of fictitious microstates of the chart of transitions on which physical conditions of QS are split. At a stage of summation of probabilities of microstates of circles their complex parts are annihilated also result of calculation - probability of number of requests in system - becomes material. Distribution of number of applications and Laplace transformation of a waiting time is received.

In work the results of calculation of QS received by means of R-approximation and by means of imitating modeling are compared. Various initial distributions of cooling - determined uniform, scale with form parameter 0,5, exponential and Weibull with form parameter 0,46 are considered. It is shown that the accuracy of calculation of stationary distribution of number of requests and the first three moments waiting time is well. The maximum distance of Kolmogorov (for Weibull distribution) made 0,043. It is shown that when coefficient of a variation of cooling time distribution is increased, the average time of expectation is increased too.

Scope of the offered calculation procedure is the analysis of functioning of organizational and technical systems with queues in which after end of the period of continuous employment there comes the temporary delay connected with the following reasons:

- maintenance of system;

- carrying out maintenance;

- restoration of the material resources spent at service;

- need of transfer to a higher body of information directorate about number of the requests served during continuous employment. Its expansion on a case of non-Markov service and the recurrent arrival, and also calculation of a queuing network which nodes is

QS with "cooling" can be the further direction of development of the offered calculation procedure.

REFERENCES

1. Khomonenko, A.D., Adadurov, S.E. & Gindin, S.I. (2013) Numerical calculations of multichannel queuing system with recurrent

input and "warm up". Izvestiya Peterburgskogo universiteta putey soobscheniya - Proceedings of Petersburg Transport University. 4(37). pp. 92-101. (In Russian).

2. Ryzhikov, Yu.I. (2013) Algoritmicheskiy podkhod k zadacham massovogo obsluzhivaniya [A use of algorithm approach in the queu-

ing theory]. St. Petersburg: Mozhaisk Military Space Academy.

3. Khomonenko, A.D., Khalil, M.M. & Gindin, S.I. (2016) [Modelling of cloud computing using multichannel queuing system with a

cooling]. Proc. of The IXI Intern. Conference on Soft Computing and Measurements. St. Petersburg. SCM'2016. St. Petersburg. May 25-27, 2016. pp. 247-251. (In Russian).

4. Kendall, M. & Stewart, A. (1966) Teoriya raspredeleniy [The advanced theory of statistics. Distribution theory]. Translated from

English. Mosow: Nauka.

5. Cox, D.R. (1955) A Use of Complex Probabilities in the Theory of Stochastic Process. Proc. of the Cambridge Phil. Soc. pp. 313-

323. DOI: 10.1017/S0305004100030231

6. Ryzhikov, Yu.I. & Ulanov, A.V. (2014) The method of calculating M/H2/n-H2 queuing system with impatient customers. Vestnik

Tomskogo gosudarstvennogo universiteta. Upravlenie, vychislitel'naya tekhnika i informatika - Tomsk State University Journal of Control and Computer Science. 2(27). pp. 47-53. (In Russian).

7. Ryzhikov, Yu.I. & Lokhvitsky, V.A. (2012) [The method of calculating of G-network]. Informatsionnye tekhnologii v upravlenii

[Information technology in management]. Proc. of the Conference. St. Petersburg. pp. 354-359. (In Russian).

8. Ryzhikov, Yu.I. & Ulanov, A.V. (2016) A use of hyperexponential distribution in non-markovian queuing systems analyses. Vestnik

Tomskogo gosudarstvennogo universiteta. Upravlenie, vychislitel'naya tekhnika i informatika - Tomsk State University Journal of Control and Computer Science. 3(36). pp. 60-65. (In Russian). DOI: 10.17223/19988605/36/6.

9. Takahashi, Y. & Takami, Y. (1976) A Numerical Method for the Steady-State Probabilities of a GI/G/c Queuing System in a General

Class. J. of the Operat. Res. Soc. of Japan. 19(2). pp. 147-157.

10. Ryzhikov, Yu.I. & Khomonenko, A.D. (1980) Iteratsionnyy metod rascheta mnogokanal'nykh sistem s proizvol'nym zakonom ob-sluzhivaniya [Iterative methods of multichannel queuing system calculation with arbitrary service time distribution]. Problemy up-ravleniya i teorii informatsii - Problems of Control and Information Theory. 9(3). pp. 203-213.

11. Serfozo, R. (2009) Basics of Applied Stochastic Processes. Berlin: Springer.