Имитационное моделирование неоднородное© процесса деерадации в системе с постепенными и внезапными отказами Текст научной статьи по специальности «Математика»

Труды Карельского научного центра РАН № 7. 2018. С. 3-13 DOI: 10.17076/mat836

УДК 519.218.3 (519.23)

ИМИТАЦИОННОЕ МОДЕЛИРОВАНИЕ НЕОДНОРОДНОГО ПРОЦЕССА ДЕГРАДАЦИИ В СИСТЕМЕ С ПОСТЕПЕННЫМИ И ВНЕЗАПНЫМИ ОТКАЗАМИ

А. В. Бородина1, В. А. Тищенко2

1 Институт прикладных математических исследований КарНЦ РАН, ФИЦ «Карельский научный центр РАН», Петрозаводск, Россия 2Петрозаводский государственный университет, Россия

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

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

A. V. Borodina, V. A. Tishenko. SIMULATION OF A NONHOMOGENEOUS DEGRADATION PROCESS IN A SYSTEM WITH GRADUAL AND INSTANTANEOUS FAILURES

The paper deals with a heterogeneous degradation process in a system with gradual and instantaneous failures. For the case when the times in the degradation stages are independent and distributed exponentially with different parameters, analytical expressions are obtained for the probability of a failure on the regeneration cycle, average time until failure on the cycle, average cycle length with and without failure, average total cycle length. For the case of arbitrary distributions in highly reliable systems, where failure is a rare event, a simulation splitting algorithm is proposed to speed up the construction of regeneration cycles. The results of experiments obtained by the splitting and the Monte Carlo methods are presented and compared with analytical solutions.

Keywords: degradation process; system with failures; probability of failure; rare events; splitting method; regenerative method; nonhomogeneous case.

0

Введение

В настоящее время изучение процессов старения и деградации в технических и биологических системах является актуальной задачей и предполагает исследование систем с постепенными отказами, для которых было разработано большое количество многоуровневых моделей надежности (см., например, [15, 18]).

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

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

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

Рассмотрим процесс деградации X = {X(¿)}^о с пространством состояний Е = {0,1,..., Ь,..., М,..., К, Е}, которые представляют стадии деградации системы (см. рис. 1).

Подробно такой процесс описан в работе [7], где для исследования динамики состояния антикоррозийного покрытия рассматривается система, управляемая двухуровневой политикой (К,Ь), и решается задача поиска оптимального управления профилактическим ремонтом с целью минимизации средних затрат при заданной системе штрафов.

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

роятность отказа системы очень мала, либо число стадий достаточно велико и т. д.

Далее будем считать, что политика управления, определяющая, когда осуществляется профилактический ремонт и когда проводить восстановление (частичное или полное), фиксирована и значения (К, Ь) являются известными.

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

Процесс стартует из состояния X(0) =0 и затем последовательно проходит К — 1 стадию деградации до перехода на стадию профилактического ремонта с номером К, если внезапный отказ не произошел. Внезапный отказ возможен начиная со стадии М, когда система не успевает перейти на стадию К, т. е. время перехода со стадии М до стадии К превышает время до отказа V, при этом V является случайной величиной (с.в.) с известным законом распределения.

Обозначим время перехода из состояния ] в состояние ] + 1 случайной величиной (с.в.) Т3 и будем считать все Tj независимыми. Если {Т} одинаково распределены (н.о.р.), то процесс называют однородным, если возможны разные распределения — неоднородным.

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

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

Определим суммарное время перехода из состояния М в состояние К значением случайной величины к_ 1

БМК = ^ Т, (1)

j=M

тогда внезапный отказ связан с событием {вык ^ V}. Предполагается, что с.в. Бык и V независимы. Если произошел отказ, то система становится неработоспособной, для полного восстановления требуется время Пр, и затем процесс заново стартует из начального состояния 0. Если отказа не случилось, то после профилактического ремонта за время 11кь система возвращается в состояние Ь (стадия возврата).

©

Рис. 1. Динамика процесса деградации с двумя типами циклов регенерации Fig. 1. Dynamics of the degradation process with two types of regeneration cycles

Регенеративная структура процесса деградации

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

Тп+1 = М> Тп : X= М},

где п ^ 0, То := 0, ^ - это момент перехода на стадию к ^ 1, а длины циклов Ук = Тк+1 — Тк, к ^ 1, являются н.о.р.

Введем необходимые обозначения для функций

Ъ(г) = Р(т < г); Гу(г) = Р(У < г); (г) = Р(Бгз < г) = * Г(г) =

Г Fitj-i(i - v)dFj(v); J 0-

Рг3 (г) = 1 — Гг3 (г)

Гг,г+1(г) = Гг(г), Гц(г) = 0,

где * означает свертку.

В работе требуется вычислить вероятность внезапного отказа на цикле регенерации

рр = Р(Бык > V) = Гмк(гушу(г), (2)

о

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

Кроме того, одним из важнейших показателей надежности является функция надежности:

я(г) = Р(т > г\х(0) = 0), г ^ о, (3)

где Т - время жизни системы (время безотказной работы). Даже в случае, когда стадии деградации распределены экспоненциально, вычисление функции надежности затруднительно, т. к. связано с нахождением обратного преобразования Лапласа (подробнее см. в [5]). Для надежных систем с редкими отказами вместо функции надежности (3) используют асимптотическую функцию надежности

PF

Ra(t) = e 1 ],

(4)

где значения рр и ] неизвестны. В таких

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

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

Ур = V + ир + Бо,м

умр = Бмк + икь + Брм, (5)

где с.в. V, UF, Som = Y1 fl/ T

j=0

независимы

лМ- 1 ,

так же, как и с.в. Бмк, икь, Брм = ^м=ь Т. Тогда безусловная (типичная) длина цикла регенерации определяется выражением

у = ур ■ 1{у^мк} + УМР ■ !{Ямк<у}. (6)

оо

Из (6) следует, что

Е[У ] = Е [тт{^ Бык }] + (Е [Пр ]+ Е[Боы ])ре

+ (Е ркь]+ Е [Бьы ])(1 — Ре ), (7)

где

E[min{V, Бык }] = / Еу (Ь)Еык (г)сН. о

Обозначим Тр время до внезапного отказа на цикле в случае, когда отказ произошел. Пусть существует Еу(ж) = /у(ж), тогда

Е[Тр] = Е^IV < Бык]

Г™ у /у (у) Р(Бык ^ У)

Уо P(V < Бык)

1 /" ™ _

= — У/у (у)Еык (y)dy. Ре 7о

Средняя длина цикла с отказом

Е[Ур] = Е[УIV < Бык]

= Е[Тр ]+ Е [Пр ]+ Е[Боы ].

йу

(9)

При условии, что существует Е3мк (ж) := /я(ж), средняя длина цикла без отказа определяется выражением

Е[У^р] = Е[У|Бык < V]

1 Г™ _

= 1- У/я (у) ЕУ (у) йу

1 — Ре 7о

+ Е[Пкь] + Е [Бьы ].

(10)

Для получения аналитических выражений характеристик, описанных выше, необходимо вычислить свертки, что в общем случае сделать практически невозможно. Более того, для нахождения явного решения необходимо выполнить инверсию преобразования Лапласа - Стилтьеса, что возможно сделать в некоторых простейших случаях. Численное вычисление обратного преобразования Лапласа показало низкую точность, в частности для большого числа стадий (см. [5]).

Аналитические формулы подсчета важнейших характеристик для однородного процесса X в случае экспоненциальных Т были получены в работе [5]. Далее получим формулы для неоднородного процесса деградации, когда все параметры распределены экспоненциально.

Неоднородный случай с экспоненциальными стадиями деградации

Рассмотрим случай, когда независимые времена Т распределены экспоненциально с

параметрами Aj (далее этот факт будем обозначать Т ~ Ежр(Аз-)) с плотностью

/ (ж) = Ajв-х'х, ж ^ 0.

Поскольку процесс деградации, как правило, ускоряется со временем, то далее будем считать, что:

Ао < ••• < Ак-1, ,3 = [0,К — 1], (11)

где параметр V определяет экспоненциальное время до отказа V ~ Eжp(v).

Утверждение 1. Если времена пребывания {Т} ~ Ежр(Аз-) являются независимыми, а время до отказа V ~ Eжp(v), то справедливо

к 1 к 1

чад = ре П лА

1

3=ы 3=ы П з

Аз (Аз + v)2Aj

(12)

где Аз = Л^ып^(Аг — ) и вероятность отказа

к-1 к-1

Ре = п а 2

з=ы з=ы

* =ы А (А + V)Аз'

(13)

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

(Т, (¿) = (1 — ^ )

1

(14)

В силу независимости Т характеристическая функция (Ямк(¿) := (ядля суммы Бык является произведением соответствующих характеристических функций

(^ = Еей(Тм + +Тк—1)

= Ее

ИТм

■ ■■Ее1гТк—1

к-1 Н

п (1—£ )-1.

(15)

По индукции можно показать (см., например, [2]), что справедливо

к1

к1

(я (*)=£ (т, (V П

Аг

з=ы

г=ы,г=з

, Аг — Аз

(16)

Сделав замену в формуле (16), к-1

П (Аг — Аз )= Аз 1=ы,г=

(17)

оо

V

получим следующее представление к 1 к 1

^ (г) = П >3 Е ^т, (г)

3=м з=м

1

>3 А''

(18)

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

1 С ж

¡з (Х) = 2П] №, (19)

из тождества (18) немедленно получим к 1 к 1

¡з (х) = П >3'Е ¡т, (х)

3=м 3=м

1

Л А'

к 1 к 1

л

3Е а

— Л, х

(20)

3=м з=м 3

тогда

__Г+Ж

Гмк(х) = Р(Бмк >х)= ¡я(у)йу =

х

к-1 к-1 Р-Л3 х

П>>Е ЛА ■

3=м 3=м 3 3

(21)

к-1 к-1

11 к 1 1 „ж

рр = - П Л3 Е 33 0 ^+">У(22)

3=м 3=м 3 3 о

и для икь ~ Ехр(ц) известно Щикь] = 1/ц, то на основании формул (10) и (20) можно сформулировать следующее утверждение.

Утверждение 2. Если времена пребывания Тз ~ Ехр(Лз) на стадиях деградации являются независимыми, время до оказа V ~ Ехр(и), а время профилактического ремонта икь ~ Ехр(ц), то справедливо

ЩУмр ] =

1

к-1 к-1

П >зЕ

1

(1 — рр) 11 "3 ^ (Л3 + и)2А3

— р 3=м 3=м 3 3

+

1 м-11

+ - + > Т-.

ц = >з

3=ь

(23)

где Аз = П^=^[1=3(>1 — >3) и вероятность от-

каза

-1 -1

рр = П >3 Е

— Л3 (Л3 + V)А3'

3=м 3=м 3 3 3

Кроме того, на основании Утверждения 1 и выражения (9) при ир ~ Ехр(цр) справедлива формула вычисления средней длины цикла с отказом:

1 м-1 1 Е[Ур ] = Е[Тр ] , (24)

/ / 771 ' * а

Поскольку с.в. V ~ Ехр(и) с плотностью ¡у (х) = ив-их, х ^ 0, то согласно формулам (2) и (21)

откуда немедленно следует выражение (13). (Альтернативный способ вывода формулы (13) приведен в статье [7], где используется свойство марковости рассматриваемого процесса и с применением преобразования Лапласа решается система дифференциальных уравнений Колмогорова.)

Аналогичным образом из формулы (8), применяя интегрирование по частям, легко получим (12).

□

3=о

где Е[Тр] вычисляется по формуле (12). Поскольку

Е[шт(^ £м к}] = Рр/V,

то, следуя приведенным выше рассуждениям, из формулы (7) получим, что средняя (безусловная) длина цикла равна

1 1 м-1 1

е[у ] = ( й + + Е т)рр +

V Цр *—' >3

р 3=о 3

1 м-1 1 + (7,+ Е >-)(1 — рр).

3=ь

(25)

Учитывая, что

м -1

Е[Бьм] = Е Лт

3=ь Л3

V

Имитационное моделирование

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

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

В статье [5] были представлены результаты экспериментов для однородного случая с экспоненциальными стадиями деградации, полученные тремя способами: аналитически, методом Монте-Карло (MC, Monte Carlo), методом расщепления (RS, Regenerative Splitting).

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

Эксперименты проводились с использованием двух языков программирования python и Си++ (в том числе с целью проверить, насколько интерпретируемый язык проигрывает в скорости получения оценки) с использованием ЭВМ: Intel(R) Xeon(R) CPU E5-2630 2.30GHz, 4GB RAM, операционная система Linux openSUSE 42.2.

Метод Монте-Карло для процесса деградации

Сначала рассмотрим стандартный метод Монте-Карло. Очевидно, что для оценивания

вероятности отказа pf на цикле регенерации (в частности это следует из формулы (13)) достаточно моделировать процесс деградации от стадии M до стадии K — 1. Таким образом, метод Монте-Карло для оценивания вероятности Pf сводится к простому алгоритму моделирования суммы независимых с.в.

Для определения других характеристик процесса необходимо имитационное построение полного цикла регенерации между моментами прихода на стадию M. В этом случае, используя регенеративный подход (см., например, [10-12]) и определение цикла (5) и (6), можно построить точечные и интервальные оценки для средней длины цикла, средней длины цикла с отказом и без, среднего времени до отказа на цикле и т. д.

Зафиксируем параметры модели v = 0, 5, ^F = 1,5, ^ = 2, L = 1, M = 5, K = 17. Для наблюдения за относительной ошибкой оценивания (relative error, [20]), вычисляемой по формуле

RE[Pf ] =

л/VarlpF ] E[Pf ] :

(26)

в экспериментах будем варьировать число циклов регенерации n и последовательность величин

X,- = Xk-i — (K — j)s, j e [0, K — 2],

где Ak-i будем выбирать, а остальные значения сдвигать на шаг s, таким образом гарантируется выполнение условия (11).

В таблице 1 приведены оценки вероятности отказа в сравнении с результатами, полученными по формуле (13). Время моделирования указано в секундах для реализации на языке python.

Таблица 1. Оценивание pF методом Монте-Карло (python) Table 1. Estimation of pF by Monte Carlo method (python)

Ak-i, s n Pf Pfmc tua Varua REua

103, 50 104 8, 75■10-3 8, 79 ■ 10-3 0,16 7,18 ■ 10-7 0, 096

104, 5 ■ 102 104 8, 79■10-4 8, 85■10-4 0,17 1,03 ■ 10-7 0, 363

105, 5 ■ 103 105 8, 79 ■ 10-5 8, 54■10-5 1, 78 1,20■10-9 0, 406

106, 5 ■ 104 107 8, 79 ■ 10-6 8, 85■10-6 205,14 1, 64 ■ 10-i2 0,145

107, 5 ■ 105 108 8, 79 ■ 10-7 8, 68 ■ 10-7 1983,11 9, 00 ■ 10-i5 0,109

Отметим, что для величин pF порядка 10-5 и меньше при числе циклов менее п = 106 метод Монте-Карло выдает pF = 0.

Замечание 1. Поскольку для метода Монте-Карло несмещенная оценка вероятности отказа pF = 1]/п, где индикатор ] = 1 для цикла с отказом, то при малых значениях

©

pF ^ 0 из (26) следует

remc =

л/рЁ^Л—рЁУ/п

PF

л/npF

В таблице 1 указана статистическая оценка КЕмс, для вычисления которой в формуле (26) использовались несмещенные статистические оценки дисперсии и математического ожидания, вычисляемые по выборке из

50 значений. Как видно из таблицы 1, общее время эксперимента для вероятности порядка 10-7 достаточно велико и составляет около 27, 5 часа.

В таблице 2 можно проследить, как меняется время вычисления оценки с уменьшением относительной ошибки оценивания (за счет увеличения числа циклов п), остальные параметры те же, что и для результатов из таблицы 1.

1

Таблица 2. Увеличение времени оценивания pp с ростом n методом Монте-Карло (python) Table 2. An increase in the estimation time of pF by the Monte Carlo method (python)

PF n pfmc tMC REmc

8, 75 ■10-3 105 8, 83 ■ 10-3 1, 71 0, 025

5 ■ 105 8, 73 ■ 10-3 9,85 0, 013

8,79 ■10-4 105 9, 00 ■ 10-4 1, 76 0,145

5 ■ 105 8, 77 ■ 10-4 8, 91 0, 040

8, 79 ■ 10-5 5 ■ 105 9, 32 ■ 10-5 8, 65 0,142

106 8, 58 ■ 10-5 19, 73 0, 085

Таким образом, с увеличением точности вычисления оценки в методе Монте-Карло существенно возрастает время моделирования. В случаях, когда аналитическое решение получить невозможно, необходимо рассмотреть алгоритмы для ускоренного имитационного моделирования циклов регенерации.

Метод расщепления для процесса деградации

В работе [5] для однородного процесса деградации в случае редких отказов был предложен имитационный алгоритм на основе динамического метода расщепления (см. [9, 13, 14]).

Ключевая проблема — это случайность порогового значения времени до отказа V, тогда как стандартная постановка задачи для ускоренного метода расщепления предполагает наличие фиксированного порога, превышение которого является редким событием.

Необходимо отметить, что и для других ускоренных методов, в частности методов уменьшения дисперсии оценки, пороговое значение в задаче оценивания вероятностей редких событий также предполагается константой (см., например, RESTART [21], условный метод Монте-Карло [6]). Для моделирования циклов регенерации в деградирующей системе техника расщепления подходит наилучшим образом и позволяет оценивать сразу несколько характеристик за один прогон. Более того, для процесса деградации возможно комбинировать расщепление с условным методом

Монте-Карло с целью уменьшения дисперсии оценки.

Другая проблема заключается в подборе оптимальных параметров: количества уровней, самих уровней и числа расщеплений на каждом уровне. Как и для однородного процесса деградации в [5], система уровней для расщепления жестко определена и совпадает со стадиями деградации. Поскольку пороговое значение V не константа, то оптимальные параметры для процедуры расщепления возможно подобрать только экспериментально.

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

При достижении процессом стадии г происходит расщепление траектории, когда для с.в. Тг генерируется Кг реализаций (копий), М ^ г ^ К — 1. Определим Км = 1 и будем считать, что с уровня М стартует одна траектория, которая разветвляется при расщеплении, и новые Кг траекторий после стадии г развиваются независимо, но имеют общую предысторию.

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

Б = Км • • • Кк-1

зависимых циклов, и мы получаем Б реализаций с.в. Бмк на каждую группу. Циклы

0

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

г-Б

Аг = £ I (з\г = !,..., Км-1,

3=(г-1)-0+1

где индикатор I^ = 1 для цикла с отказом (ЮО = о, иначе) и с.в. Аг н.о.р.

Последовательность {I^ 1} является дискретной регенерирующей с постоянной длиной цикла вг = О и моментами регенерации г ■ О, г е [1,Км-1 ]. Тогда несмещенная оценка рр вероятности рр является строго состоятельной, т. е. при Км-1 ^ те с вероятностью 1 имеет место сходимость

PF =

ERm-1 л 7=1 Aj

rm- 1 ■ d

E Ej=i I(7)

D

= PF. (27)

Поскольку для экспоненциального неоднородного случая 'Pf можно вычислить явно по формуле (13), то для удобства сравнения оценки с точным решением подсчитаем относительную экспериментальную ошибку (relative experimental error, [20]) для обоих методов MC и RS по формуле (в %):

RER[pF] = |pF - PF | ■ 100/pf.

В таблице 3 приведены результаты запуска программы на языке Си++ для метода Монте-Карло и метода расщепления. Как и следовало ожидать, время моделирования методом Монте-Карло на языке СИ++ при одинаковых фиксированных параметрах существенно меньше, чем для языка python (см. табл. 1).

Таблица 3. Сравнение оценок pF в методах Монте-Карло и расщепления (Си++) Table 3. Comparison of the estimates of pF in Monte Carlo and splitting methods (C++)

Ак-i, s n PF RERmc RERrs tMC tRS REMC rers

103, 50 104 8, 75■10-3 0,42 0, 48 0,037 0, 012 0,104 0, 984

104, 5 ■ 102 104 8,79■10-4 2, 89 1, 69 0,037 0, 010 0, 318 2,167

105, 5 ■ 103 105 8, 79 ■ 10-5 3, 84 0, 95 0, 389 0, 095 0, 366 3, 254

106, 5 ■ 104 107 8, 79 ■ 10-6 0, 75 1, 80 37, 83 6, 56 0,125 1, 967

107, 5 ■ 105 108 8, 79 ■ 10-7 8, 51 3, 99 378, 0 25, 06 0,166 3,188

Эксперименты показали, что, как и в однородном случае, метод расщепления дает выигрыш по времени построения циклов регенерации. Более того, из таблицы 3 видно (столбцы КЕКмс, КЕКяя), что оценка по методу расщепления является более близкой к аналитическому решению, т. к. имеет меньшую относительную экспериментальную ошибку. Более того, для метода расщепления во всех экспериментах КЕК < 4 %, тогда как метод Монте-Карло показал КЕК > 8 % для вероятности порядка 10-7. Тем не менее для заданного числа циклов метод расщепления имеет большую дисперсию, чем метод Монте-Карло, а следовательно, и большую относительную ошибку оценивания.

Замечание 2. Качество оценивания в методе расщепления существенно зависит от выбора уровней и числа расщеплений. Тем не ме-

нее для процесса деградации система уровней определена изначально стадиями деградации процесса и не может подбираться оптимально. Поскольку время до отказа V является случайной величиной, то использовать процедуру настройки как для стандартного метода расщепления (pilot run [17, 19]) не представляется возможным. В этой связи для уменьшения относительной ошибки оценивания возможно комбинировать метод расщепления с условным методом Монте-Карло [6].

Для вычисления асимптотической функции надежности (4) необходимо получить оценку средней длины цикла без отказа E[YNf]. Результаты моделирования методами Монте-Карло и расщепления в таблице 4 можно сравнить с аналитическим выражением по формуле (23).

Таблица 4- Результаты вычисления E[YNf] методами MC и RS Table 4- Estimation results of E[YNF] by MC and RS methods

ak-i n e[ynf ] e[ynf ]mc e[ynf ]rs

103 105 0,528497 0,528150 0,528275

104 105 0,502851 0,502976 0,502807

105 105 0,500285 0,500063 0,500282

106 107 0,500029 0,499987 0,500030

107 108 0,500003 0,499961 0,500438

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

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

Принцип построения доверительного интервала для процесса деградации с использованием метода расщепления был описан в работе [5]. Напомним, что классический регенеративный метод предполагает независимость и одинаковую распределенность циклов регенерации.

Когда для построения циклов используется метод расщепления, то очевидным образом циклы регенерации являются зависимыми (подробнее см. в [4]).

Для случая зависимых, но одинаково распределенных циклов можно воспользоваться центральной предельной теоремой для стационарной в узком смысле последовательности к-зависимых с.в. с конечной дисперсией (см. [1, 8]). При к = 1 для однозависимых циклов формула доверительного оценивания содержит ковариацию между двумя зависимыми циклами (см. в [16]).

Заключение

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

Финансовое обеспечение исследований осуществлялось из средств федерального бюджета на выполнение государственного задания КарНЦ РАН (Институт прикладных математических исследований КарНЦ РАН) и при финансовой поддержке РФФИ (проекты 18-07-00187, 18-07-00147).

Литература

1. Ширяев А. Н. Вероятность. М.: Наука, 1980. 576 с.

2. Bibinger M. Notes on the sum and maximum of independent exponentially distributed random variables with different scale parameters. arXiv preprint arXiv:1307.3945, 2013.

3. Botev Z. I, Kroese D. P. Efficient Monte Carlo simulation via the generalized splitting method // Statistics and Computing. 2012. Vol. 22, no. 1. P. 1-16. doi: 10.1007/s11222-010-9201-4

4. Borodina A. V. Regenerative modification of the splitting method for estimating the overload probability in queuing systems: PhD Thesis. Petrozavodsk State University, 2008 (in Russian).

5. Borodina A. VEfrosinin D. VMorozov E. V. Application of Splitting to Failure Estimation in Controllable Degradation System // Vishnevskiy V., Samouylov K., Kozyrev D. (eds). Distributed Computer and Communication Networks. Communications in Computer and Information Science, Springer, Cham. 2017. Vol. 700. P. 217-230. doi: 10.1007/978-3-319-66836-9

6. Chan J. C. CKroese D. P. Rare-event probability estimation with conditional Monte Carlo // Ann. Oper. Res. 2011. Vol. 189. P. 43. doi: 10.1007/s10479-009-0539-y

7. Efrosinin D. V., Farhadov M. P. Optimal management of the system with the gradual and instantaneous failures // Dependability. 2009. No. 1(28). P. 27-42 (in Russian).

8. Ferguson T. S. A course in large sample theory: Texts in Statistical Science. Chapman and Hall/CRC, 1996. 256 p.

9. Carvels M. The splitting method in rare event simulation: PhD Thesis. The University of Twente, The Netherlands, 2000.

10. Glynn P. W. Some topics in regenerative steady state Simulation // Acta Applic. Math. 1994. Vol. 34. P. 225-236. doi: 10.1007/BF00994267

11. Glynn P. W, Iglehart D. L. Conditions for the applicability of the regenerative method // Management Science. 1993. Vol. 39. P. 1108-1111. doi: 10.1287/mnsc.39.9.1108

12. Glynn P. W., Iglehart D. L. A joint central limit theorem for the sample mean and regenerative variance estimator // Annals of Operations Research. 1987. Vol. 8. P. 41-55. doi: 10.1007/BF02187081

13. Heegaard P. E. A survey of Speedup simulation Techniques: Workshop tutorial on Rare Event Simulation, Aachen, Germany, 1997.

14. Heidelberger P. Fast simulation of rare events in queuieng and relaibility Models // Performance Evaluation of Computers and Communications Systems. Springer-Verlag, LN in Computer Sci, 1993. Vol. 729. P. 165-202.

15. Lisnuansky A., Levitin G. Multi-state system reliability: assessment, optimization and application. New Jersey, London, Singapore, HongKong: World Scientific, 2003. doi: 10.1142/5221

16. Morozov E., Aminova I. Steady-state simulation of some weak regenerative networks

// European Transactions on Telecommunications ETT. 2002. Vol. 13, no. 4. P. 409-418. doi: 10.1002/ett.4460130412

17. Rubinstein R. Y, Kroese D. P. Simulation and the Monte Carlo method. Hoboken: John Wiley & Sons, Inc., 2016. 396 p. doi: 10.1002/9781118631980

18. Rykov V., Dimitrov B. On multi-state reliability systems // Proc. of Seminar Applied Stochastic Models and Information Processes. 2002. P. 128-135.

19. Rubinstein R. Y., Ridder A., Vaisman R. Fast Sequential Monte Carlo Methods for Counting and Optimization. Hoboken: John Wiley & Sons, Inc., 2014. 208 p. doi: 10.1002/9781118612323

20. Vaisman R., Roughan M., Kroese D. P. The multilevel splitting algorithm for graph coloring with application to the Potts model // Philosophical Magazine. 2017. doi: 10.1080/14786435.2017.1312023

21. Vill'en Altamirano J. RESTART Vs Splitting: A Comparative Study // Proceedings of the 11th Workshop on Rare Event Simulation, RESIM'16, Eindhoven (The Netherlands), 2017. P. 1-12. doi: 10.1016/j.peva.2018.02.002

Поступила в редакцию 08.04.2018

References

1. Shiryaev A. N. Veroyatnost' [Probability]. Moscow: Nauka, 1980.

2. Bibinger M. Notes on the sum and maximum of independent exponentially distributed random variables with different scale parameters. arXiv preprint arXiv:1307.3945, 2013.

3. Botev Z. I., Kroese D. P. Efficient Monte Carlo simulation via the generalized splitting method. Statistics and Computing. 2012. Vol. 22, no. 1. P. 1-16. doi: 10.1007/s11222-010-9201-4

4. Borodina A. V. Regenerative modification of the splitting method for estimating the overload probability in queuing systems: PhD Thesis. Petrozavodsk State University, 2008 (in Russian).

5. Borodina A. V., Efrosinin D. V., Morozov E. V. Application of Splitting to Failure Estimation in Controllable Degradation System. In: Vishnevskiy V., Samouylov K., Kozyrev D. (eds). Distributed Computer and Communication Networks. Communications in Computer and Information Science, Springer, Cham. 2017. Vol. 700. P. 217-230. doi: 10.1007/978-3-319-66836-9

6. Chan J. C. C, Kroese D. P. Rare-event probability estimation with conditional Monte Carlo. Ann. Oper. Res. 2011. Vol. 189. P. 43. doi: 10.1007/s10479-009-0539-y

7. Efrosinin D. V., Farhadov M. P. Optimal management of the system with the gradual and instantaneous failures. Dependability. 2009. No. 1(28). P. 27-42 (in Russian).

8. Ferguson T. S. A course in large sample theory. Texts in Statistical Science. Chapman and Hall/CRC, 1996. 256 p.

9. Garvels M. The splitting method in rare event simulation: PhD Thesis. The University of Twente, The Netherlands, 2000.

10. Glynn P. W. Some topics in regenerative steady state simulation. Acta Applic. Math. 1994. Vol. 34. P. 225-236. doi: 10.1007/BF00994267

11. Glynn P. W, Iglehart D. L. Conditions for the applicability of the regenerative method. Management Science. 1993. Vol. 39. P. 1108-1111. doi: 10.1287/mnsc.39.9.1108

12. Glynn P. W, Iglehart D. L. A joint central limit theorem for the sample mean and regenerative variance estimator. Annals of Operations Research. 1987. Vol. 8. P. 41-55. doi: 10.1007/BF02187081

13. Heegaard P. E. A survey of Speedup simulation techniques. Workshop tutorial on Rare Event Simulation, Aachen, Germany, 1997.

14. Heidelberger P. Fast simulation of rare events in queuieng and relaibility Models. Performance Evaluation of Computers and Communications

Systems. Springer-Verlag, LN in Computer Sci., 1993. Vol. 729. P. 165-202.

Stochastic Models and Information Processes. 2002. P. 128-135.

15. Lisnuansky A., Levitin G. Multi-state system reliability: assessment, optimization and application. New Jersey, London, Singapore, HongKong: World Scientific, 2003. doi: 10.1142/5221

16. Morozov E., Aminova I. Steady-state simulation of some weak regenerative networks. European Transactions on Telecommunications. 2002. Vol. 13, no. 4. P. 409-418. doi: 10.1002/ett.4460130412

17. Rubinstein R. Y, Kroese D. P. Simulation and the Monte Carlo method. Hoboken: John Wiley & Sons, Inc., 2016. 396 p. doi: 10.1002/9781118631980

18. Rykov V., Dimitrov B. On multi-state reliability systems. Proc. of Seminar Applied

19. Rubinstein R. Y., Ridder A., Vaisman R . Fast Sequential Monte Carlo Methods for Counting and Optimization. Hoboken: John Wiley & Sons, Inc., 2014. 208 p. doi: 10.1002/9781118612323

20. Vaisman R., Roughan M., Kroese D. P. The multilevel splitting algorithm for graph coloring with application to the Potts model. Philosophical Magazine. 2017. doi: 10.1080/14786435.2017.1312023

21. Vill'en Altamirano J. RESTART Vs Splitting: A Comparative Study. Proceedings of the 11th Workshop on Rare Event Simulation, RESIM'16, Eindhoven (The Netherlands), 2017. P. 1-12. doi: 10.1016/j.peva.2018.02.002

Received April 8, 2018

СВЕДЕНИЯ ОБ АВТОРАХ:

Бородина Александра Валентиновна

научный сотрудник, к. ф.-м. н.

Институт прикладных математических исследований

КарНЦ РАН, Федеральный исследовательский центр

«Карельский научный центр РАН»

ул. Пушкинская, 11, Петрозаводск,

Республика Карелия, Россия, 185910

эл. почта: borodina@krc.karelia.ru

тел.: (8142) 763370

Тищенко Виталий Анатольевич

студент

Петрозаводский государственный университет

пр. Ленина, 33, Петрозаводск, Республика Карелия,

Россия, 185910

эл. почта: vitalik1tishenko@gmail.com тел.: (8142) 719606

CONTRIBUTORS:

Borodina, Alexandra

Institute of Applied Mathematical Research, Karelian Research Centre, Russian Academy of Science 11 Pushkinskaya St., 185910 Petrozavodsk, Karelia, Russia e-mail: borodina@krc.karelia.ru tel.: (8142) 763370

Tishenko, Vitalii

Petrozavodsk State University

33 Lenin St., 185910 Petrozavodsk, Karelia, Russia

e-mail: vitalik1tishenko@gmail.com

tel.: (8142) 719606