Аппроксимации функции восстановления и стратегия управления эксплуатационными затратами Текст научной статьи по специальности «Математика»

а-

УДК 519.218.4:62.5

DOI:https:/^doi.org/10.25728/pu.2018.4.4

АППРОКСИМАЦИИ ФУНКЦИИ ВОССТАНОВЛЕНИЯ И СТРАТЕГИЯ УПРАВЛЕНИЯ ЭКСПЛУАТАЦИОННЫМИ ЗАТРАТАМИ

В.Н. Русев, А.В. Скориков

Рассмотрена аппроксимация функции восстановления для рекуррентных потоков восстановлений. Предполагается, что время безотказной работы объектов системы описывается с помощью двухпараметрического распределения Вейбулла — Гнеденко. Предложены аналитический и численный подходы (на основе методов дискретизации) решения уравнения Вольтерра. Полученные соотношения и алгоритмы проверены на контрольных примерах с помощью профессионального математического пакета Wolfram Mathe-matica. Изложенные результаты применены в стратегии «блоковой (групповой) политики замен» для обеспечения надежного и эффективного функционирования технологических объектов.

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

ВВЕДЕНИЕ

Функция восстановления Н(?) определяется как ожидаемое число отказов системы (или элемента) в течение времени (0, ?). Этот параметр используется для построения системы менеджмента надежности, определения оптимального плана профилактического обслуживания или замен [1—3]. Функция восстановления имеет важное значение для разработки политики гарантий [4], а также применяется в принятии решений планирования схемы поставок [5], построения планов страхования [6].

Для рекуррентных потоков восстановлений интегральное уравнение Вольтерра устанавливает соотношение между функцией восстановления Н(?) и функцией распределения F(t) времени работы между отказами [7], когда плотность распределения Д?) = /"(?) определена:

H(t) = F(t) + \H(x)f(t - x)dx.

(1)

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

известна [7] асимптотическая формула для функции восстановления (асимптота Смита)

H(t) ~ 11 + 2 ц 2

2

( 2 : а - ц

2

V ц

t ^ да,

(2)

где ц и ст — математическое ожидание и дисперсия времени работы между отказами.

Отметим также публикации [8, 9], посвященные прямому вычислению функции надежности для систем в терминах сверток.

В настоящей статье предполагается, что весь жизненный цикл функционирования объектов системы описывается с помощью двухпараметричес-кого закона распределения Вейбулла — Гнеденко (см. например, работу [10]). Распределение Вейбул-ла — Гнеденко занимает важное место среди распределений времени безотказной работы систем, состоящих из групп большого числа элементов, отказы которых происходят взаимно независимо, так что отказ любого из элементов приводит к отказу всей системы (принцип «наислабейшего звена»).

Согласно государственному стандарту [11], «Планово-профилактические ремонты или замены полезны в случаях, когда отказы одной или нескольких ключевых составных частей изделия имеют четко выраженный износовый и/или усталостный характер, что соответствует описанию

о

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

Функция распределения и плотность распределения Вейбулла — Гнеденко имеют, соответственно, вид:

Щ = 11 - е-(<")Р,г > 0, 10, г < 0

f(t) = / в *-У о, t < о

1 -(atf

e ( '

t > 0

с параметрами формы в > 0, масштаба а > 0.

Авторы работы [12] W.L. Smith and M.R. Lead-better предложили метод для вычисления функции восстановления для распределения Вейбулла — Гнеденко с помощью разложения в степенной ряд для в > 1, численное вычисление этого ряда ограничено малым диапазоном t; A.G. Constantine, N.I. Robinson [13] представили метод вычисления функции восстановления H(t) с помощью вычетов, позволяющих получить представление для нее равномерно сходящимся рядом затухающих экспоненциальных членов. В статье [14] была предложена оценка функции восстановления, основанная на первых трех моментах распределения, уточняющая асимптотическую формулу (2). Получен также ряд других приближений, таких как полученных авторами работ [15—17], изучавших функцию восстановления для нормального, гамма-, равномерного распределений времени жизни.

1. АНАЛИТИЧЕСКОЕ ПРЕДСТАВЛЕНИЕ ФУНКЦИИ ВОССТАНОВЛЕНИЯ.

МЕТОД ПРОИЗВОДЯЩИХ ФУНКЦИЙ МОМЕНТОВ

Применим методы операционного исчисления. Из уравнения (1) получаем

H(s) = F (s) + H(s) f (s), откуда следует

H (s) = -Ш- = f(s) .

1 - f( s) s( 1 - f( s))

Для нахождения оригинала преобразования

Лапласа Н (л) воспользуемся производящей функции моментов плотности распределения Вейбул-ла — Гнеденко.

Соотношение к(г) = Н(г) в терминах преобразования Лапласа можно записать как Н (л) = к (л)/л.

В работах [18, 19] получено представление преобразования Лапласа плотности восстановления

h (s) = — + v1s

- ^ + 3 v 2 - 2v1 v 3 s +

2v2 J 12v1

1

+ 3 v 2 - 4 v 1 v 2 v 3 + v 1 v 4 s2 +

24v1

Отсюда

H (s) = ^ = J_ +

s v1 s2

3

r -1\ 1 + 3 v 2 - 2 v1 v 3 + 2v1 Js 12v3

1

+ 3 v 2 - 4 v 1 v 2 v 3 + v 1 v 4 s2 +

24v1

Используя обозначения

1

v1

co = -■

c = -V2

2 v1

с2 =

3 v22 - 2v1 v3

12 v 31

c, =

_ 3v3 - 4V1V2V3 + V2V4

24v41

ck =

1 0 0 0 -m0

-m0 1 0 • 0 m1

m1 -m0 1 • 0 -m2

-m1 m1 -m0 • 0 m3

m3 -m2 m1 0 1 -m4

k +1 mk 2 (-1 )kmk 3 (-1)k 1 mk 1 • • -m0 (- 1)km

где vk — начальный момент случайной величины \ с плотностью /(г) порядка к:

Vk=г ^=1+к), п а ^

r(s) = J Xs e — гамма-функция Эйлера, s е С, + (с1 — 1)1 + с2 + c3s + c4s2 + ... .

имеем

H (s) =

s

Тогда

H(t) = F(t) + L-

Ц f (s) + (q - 1)1 f (s) +

s

s

+ c2 f (s) + c3sf (s) + c/ f (s) + ...

где Ь обозначает обратное преобразование

Лапласа для функции 6(5). Таким образом,

H(t) = F(t) + L

-1

c0 f (s)

s

+ (c1 - 1) Jf(x)dx +

+ c2f(t) + f) + cf^ + ... =

= F(t) + L

-1

f (s)

+ (c1 - 1)F(t) + c2F'(t) +

+ c3F"(t) + c4F'"(t) + В силу того, что

L

-1

f (s)

= c0 J(t - x)aV - 1 e-(aT) dx =

= c0aep

t Jxe - 1 e-(aT) dx - Jxe e-(aT) dx

= c0aep

t-L y(1, (at)e) - -¿L- yV 1 + 1, (at)

aep

= С

в + 1o

a р

tF(t) - 1 yV 1 + J1' (at)1

a V в

где у(а, х) = |?а — неполная гамма-функция

Эйлера и ввиду равенств

р е-(ат) dx = 1

+1,

a р V в

(at)e

Yl^ ,(a t )в

у(1, (а?)в) = | е zdz = Д?), о

окончательно получаем представление:

tF(t) - 1 y( 1 + 1' (at)в

Н(?) = Н(?) + Со ^ ,

0 а ^ в

+ (сх - 1)Д?) + с2/"(?) + с3/"'(?) + с4/""(?) + ... = = Д?)[1 + Со?] - С0 у( 1 + 1, (а*)р) + (сх - 1)Д?) +

+ С2/"(?) + С3/"'(?) + С4/""(?) + ... =

= с°?Д?) - 0° у( 1 + в, (а0р) + с^(?) + С2^"(?) + + с3/""(?) + с4/""(?) + ..., где коэффициенты ск определены выше.

Другая запись полученного соотношения: Н(?) = Д?)(с°? + сх) - С° у( 1 + 1, (а*)Р) + С2Л<) +

+ С3/"'(?) + С4/""(?) + ...

Осталось только еще удостовериться в том, что полученное соотношение согласуется с асимптотической формой функции восстановления (2), которая в наших обозначениях имеет вид (см. далее рис. 4):

Н(?) = с°? + сх - 1 + о(1), ? ^ да.

2. ДИСКРЕТИЗАЦИЯ ИНТЕГРАЛЬНОГО УРАВНЕНИЯ

Для приближенного метода решения уравнения

(1) относительно Н(?) применим методы дискретизации интегрального уравнения, которые в вариационных задачах называются прямыми методами. К таким методам относятся конечно-разностный метод Эйлера, метод Ритца, а также известный в задачах механики метод конечных элементов.

Дискретизация по методу Галеркина с применением кубических сплайнов была рассмотрена в работе [20]. Однако счет для большого числа точек разбиений, по предложенной в этой работе схеме, затруднителен из-за большого массива данных в оперативной памяти компьютера. Численное решение уравнения восстановления (прямой Я^-ме-тод) рассмотрено в статьях [21—23], где получены алгоритмы решения уравнения методами типа метода трапеций и метода Симпсона.

Используемая далее в работе дискретизация состоит в применении соответствующего алгоритма только к неизвестной функции, стоящей под знаком интеграла, а интегралы от ядер операторов вычисляются в пакете МаШешаИса. Таким образом, решение ищется, как в методе Ритца — Галеркина, в виде линейной комбинации соответствующих базисных функций, но для определения коэффициентов, как в методе коллокаций Канторовича [24], невязки в узлах полагаются равными нулю. Точность аппроксимации решения достигается не благодаря более сложным базисным функциям, как в методе Галеркина, а благодаря более мелкому разбиению области.

Рассмотрим три способа дискретизации уравнения (1). Разделим рассматриваемый отрезок времени [0, ?тах] на п равных частей точками = 0, ?1 = + А, ..., ?п = + пА, где А — шаг разбиения.

Метод правых узлов. Будем искать решение Н(?) уравнения (1) в виде суммы кусочно-постоянных функций, равных на каждом отрезке разбиения _ 1 < ? < значению искомого решения в правой

0

c

2

s

c

s

0

0

0

0

Рис. 1. Метод правых узлов

Рис. 2. Метод средних

точке гк разбиения Нк = Н(гк), также как в квадратурном методе прямоугольников (рис. 1). Если ввести координатные функции 1к(г), равные единице при гк _ 1 < г < гк и нулю в остальных точках, то приближенное решение записывается в виде

п

Н(г) - I НкЛ(г).

к = 1

Метод средних. Положим значение приближенного решения и(г) на к-м отрезке, равным

среднему значению Нк =

_Н( гк -1) + Н( гк)

2

как в ме-

тоде трапеций квадратурных формул (рис. 2). Приближенное решение, соответственно, имеет вид:

^ п _

Н(г) = I НДО.

к = 1

Линейные сплайны. В качестве приближенного решения на к-м отрезке Нк (г) выберем интерполяционные многочлены Лагранжа первой степени, которые в узлах гк _ 1, гк, к = 1, ..., п, равны значениям решения Нк _ 1 = Н(гк _ 1), Нк = Н(гк). Таким образом, на каждом отрезке разбиения искомое решение заменяется линейной функцией, т. е. рассматривается приближение решения ломаной (рис. 3).

Найдем решение по первому способу. Подставим точки разбиения в уравнение (1), заменим Н(г) на Н(гк) = Нк. Получим:

Но = ДО) = л, <1

Н1 = ^(г1) + Н1 ¡/ц - ^ = |г1 - л = г| = о

<1

= дг^ + Н1 /»¿г, о

Н1 = + Н1КГ

<1

где = Д^), К1 = /»¿г,

о

<1 <2 Н2 = дг2) + Н1 ¡/(г2 - ^ + Н2 |дг2 - ^ = о <1

<2 <2 - <1 = |г2 - л = г| = Дг2) + Н1 | ДгМг + Н2 | Дг)^ =

<2 - <1

= дгз) + Н1 ¡/»¿г + Н2 /»¿г,

Н2 = ^2 + Н1К2 + Н2К1,

где = Дг2), К2 = ¡Дг)^,

Нп = Дг^ + Н ¡/^ - ^ + Н2 ¡/(гп - ^ + ... +

о <1

<п

+ Нп 1 /(гп - ^ = |гп - 5 = Г| = <

<п 1

= Дг^ + Н1 1 Дг)^Г + Н2 1 Дг)аТ +...+ Нп 1/>)аГ,

<п 1

<п 2

Н = Д, + #,КИ + Н,К , + ... + дх,

п п 1 п 2 п — 1 п 1'

где Рп = Дгп), Кп = 1 Дг)аТ

<п 1

Рис. 3. Метод линейных сплайнов

о

о

п

о

п

Таким образом, для п точек разбиения получаем систему линейных уравнений вида:

H0 = F0,

H = f + ад,

H2 = F2 + H1K2 + H2 Kb

H = Fn + H^ + ВД -1 + ... + ад.

(3)

n K1.

Решение системы можно записать рекуррентными формулами:

ио = Fо,

Н1 = Д / (1 - К)

и = (д + и К ) / (1 - К),

Нп = ( д + иКп — И2Кп- 1 +... + Ип- 1К2)/(1 - К1).

При втором способе дискретизации заменим Н(?) на соответствующем отрезке средним йк =

= ——. Соответственно, модифицируется система (3):

H0 = Fо,

H1 = F1+

H1+H0 -----2---------

H + H

K,

H2 = F2 + K2 + H2_+_H1 к,

1,

xx xx Hi + H0 H2 + Hi

Hn = Hn + -J—-0 Kn + -Ц—" Kn -1 +... +

+ Hn + Hn - 1 K1.

2

Отсюда получаем рекуррентные формулы для решения:

ио = Fо,

H. = V F. + F0 K.) / V1 - K

H2 = (F2 + H- K + H- K + HJo K1) / (1 - f

Hn = (Fn + H-K + f*n -1 +... + + H-*, +

+ f«. - ■ +... + ^f1 K1) / (1 - f

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

и/?) = На _ 1 Л=±. + Ик

k - 'k - 1 tk - tk - 1

Обозначим

1k(t) =

_ tk - t _ tk - t

?к - ?к - 1 А Тогда

ик(?) = ик - АО - иА - 1(?).

Для линейных функций 1к справедливо при т > к равенство

1k(tm - t) =

_ tk - tm + t _ -(tm - tk) + t _

tm - k - t —

A A

= -1m - k(t).

(4)

Из уравнения (1) ио = /(0) = далее при ? = ?1, учитывая равенство (4), имеем

'1

и = ад + |я?1 - 5)(и,/1(5) - ида)^ =

о

'1

= |?1 - 5 = г | = Д?1) + |/(г)(ио/1(?1 - Г) -

о

'1

- и^ - г)Лг = Д - ио /»^г +

о

'1

+ и1 /ЭДгМг = Д - ио^ю + и^ц. о

Отсюда

и1(1 - 611) = /1 - ио61о. Подставим в уравнение (1) ? = ?п:

'1

ип = ДО + !/(?п - *)№(*) - и^)^ +

о

'2

+ |/(?п - ^(и^) - и^/1 (5))+ ... +

'1

'я

+ 1 Д?п - 5)(и,_1/п(5) - иА-Д5)^ = |?1 - 5 = Г| =

'я 1

'я

= /1 + 1 ШЩ^п - Г) - и^ - +

'я 1

'я 1

+ 1 Аг-хи^п - Г) - - + ... +

'я 2 '1

+ 1/(г)(ип - ¿('п - Г) - ип^п - 1(?п - =

0

= F„ - Ho J f>)/„ - 1(r)dr + H1 J f(r)/„(r))dr -

^ 1 tn 1

^ 1 tn 1

- H1 J f(r)/„ - 2(r)dr + H2 J f>)/„ - 1(r)dr - ... -

tn 2

tn 2

- Hn - 1 f>)(/o(r)dr + H„ j/(r)/1(r)dr =

o

= F - HG , + H, (G - G , .) + ... +

1 0 nn - 1 1v nn n - 1n - 2'

+ Hn - 1(G22 - G10) + HnG11,

где

Gj = J f(r)/j(r)dr.

tk 1

Таким образом,

Hn(1 - G11) = F1 - H0Gn

0 nn — 1

+

+ H1(Gnn - Gn - 1n - 2) + ... + Hn - 1(G22 - G10). Итак, система уравнений имеет вид

H0 = ^ H1(1 - G11) = F1 - H0G10, H2(1 - G11) = F2 - H0 G21 + H1( G22 - G10 ),

(5)

Hn( 1 - G11) = Fn - H0Gnn - 1 + + H1( Gnn - Gn- 1n-2) + ... + Hn- 1( G22 - G10).

Рекуррентные формулы для решения системы (5) имеют вид:

u0 = f0,

H1 = (F1 - H0 G10) / (1 - G1),

H2 = С/2 - H0G21 + H1( G22 - G10))/(1 - G11),

Hn = (Fn - H0Gnn- 1 + H1( Gnn - Gn- 1n-2) + ... + + Hn - 1( G22 - G10))/(1 - Gn).

3. КОМПЬЮТЕРНОЕ МОДЕЛИРОВАНИЕ И СРАВНЕНИЕ РЕЗУЛЬТАТОВ

Полученные алгоритмы были проверены на контрольных примерах с помощью профессионального математического пакета Wolfram Mathema-tica. Погрешность приближенного решения определяется не только построенным алгоритмом, но также погрешностью численного интегрирования

в пакете Wolfram Mathematica, равной 10 15. Отметим, что погрешность решения контролируется применением различных методов и известным асимптотическим приближением (2). Для распре-

Рис. 4. Аппроксимация функции восстановления для распределения Рэлея

Рис. 5. Сравнение аппроксимаций функций восстановления для распределения Рэлея, полученных методом производящих функций моментов и методом вычетов Constantine, Robinson

деления Рэлея, важного частного случая распределения Вейбулла — Гнеденко (в = 2, а = 1), функции восстановления представлены на рис. 4. Из графиков видно, что приближение линейными сплайнами и метод средних дает одно и то же решение при малом числе точек разбиения (20 точек), которое близко к приближенному решению, полученному методом производящих функций моментов. Рассмотренные приближенные решения быстро сходятся к асимптотическому значению функции восстановления (2). Метод правых узлов для рассмотренного числа точек разбиения (20 точек) дает большую погрешность при г > 1. Однако, увеличение числа точек разбиения до 200 дает практически совпадение приближенных решений, найденных рассмотренными выше методами правых узлов, средних и производящих функций. Время счета при таком разбиении составляет приблизительно 2—3 с.

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

0

4. ОПТИМАЛЬНАЯ СТРАТЕГИЯ ПРОФИЛАКТИЧЕСКОГО ОБСЛУЖИВАНИЯ

По ГОСТ Р 27.606-2013 [11]: «Если же отказ не несет угроз безопасности, но ведет к утрате изделием готовности к применению по назначению, то периодичность замен устанавливают, исходя из заданного уровня готовности, обеспечиваемого при оптимальных затратах, включающих в себя, в том числе, стоимость заменяемых изделий и экономический ущерб от отказов». Таким образом, необходим баланс между суммой, потраченной на профилактическое обслуживание, и суммой на замены при внезапном отказе.

Применим стратегию технического обслуживания, известную как «групповую или блоковую политику замен» («block replacement policy — BRP») [1, 3]. Предполагается, что объект заменяется новым изделием при постоянной длине интервала замен независимо от возраста объекта, а также замены объекта происходят столько раз, сколько требуется на интервале (0, tp) при внезапных отказах объекта. В работе [25] рассмотрены стратегии замен по эксплуатационным затратам, в частности, уточняющие определение оптимального срока

Рис. 6. Графики приближений функции стоимости для распределения Рэлея при коэффициенте затрат с= 2,6

Рис. 7. Графики приближений функции для распределения Рэлея при коэффициенте затрат со = 10

службы при восстановлении работоспособности объекта без восполнения его ресурса.

Пусть Ср — средняя стоимость профилактического обслуживания, С/ — средняя стоимость восстановления при отказе (Ср < Су). Средняя стоимость на интервале (0, гр) профилактического обслуживания и восстановления после отказа равна Ср + СуН(гр), где Н(гр) — среднее число восстановлений на интервале (0, гр). В качестве критерия оптимальности рассматривается средняя стоимость эксплуатационных затрат (удельная стоимость затрат) единицу времени:

) = Cp( 1 + CoH tp) )

(6)

где со = Су/Ср — коэффициент затрат. Точка минимума функции (6) дает соответствующее значение времени оптимального профилактического обслуживания гр. Известно, что в точке го экстремума Л(го) = С/Я'(го) = Сук(го), т. е. значение функции Л(г) определяется плотностью восстановления к(го). Более того, можно доказать равенство

Л'Ц) = го С/#"(го). (7)

Из равенства (7) следует, что характер выпуклости функции Л(г) в точках экстремума совпадает с характером выпуклости функции Н(г) и осцилляция функции восстановления согласуется с осцилляцией средней стоимости эксплуатационных затрат.

Достаточным условием существования минимума функции Л(г) является условие на коэффициент затрат с0 [26, 27], имеющее в наших обозначениях вид:

2

Co >

1 - CV

-2 '

(8)

2 2 2

где СК = ст /ц — квадрат коэффициента вариации. Отметим, что условие (8) не является необходимым. Приведем иллюстрацию данного утверждения, служащую контрпримером условию (8). Рассмотрим распределение Рэлея, для которого ограничение (8) означает, что с0 > 2,8. Однако, как видно из графика (рис. 6) при значении коэффициента со = 2,6, не удовлетворяющему условию (8), функция Л(г) имеет вполне конкретный минимум. Более того, при этом же коэффициенте затрат функция стоимости имеет еще и максимум. Иными словами, при большем значении интервала профилактического обслуживания можно попасть в точку максимума затрат, а из формулы (2) следует, что Л(<») = С//ц является наименьшим значением функции стоимости. Поэтому в этом случае оптимальной стратегией будет замена лишь после наступления отказа.

p

Заметим, что применение критерия (6) для определения времени оптимального профилактического обслуживания требует достаточно точного вычисления функции восстановления. Сравнение графиков (рис. 7) критерия (6) показывает значительное расхождение графиков Л(?) при вычислении функции восстановления по методу производящих моментов или дискретизацией с графиком Л(?) по методу [14] с использованием трех моментов.

Соответствующие координаты точек минимума: {0,32; 6,28}; {0,32; 6,22}: {0,22; 5,17}. Сравнение результатов, полученных в данной работе с результатами {0,34; 6,17} работы [16] показывают небольшое расхождение 5 и 2 % соответственно.

Таким образом, стратегия «групповой политики восстановлений» дает возможность обеспечения надежного и эффективного функционирования технологических объектов в рамках системного подхода к мониторингу эксплуатационных показателей.

ЗАКЛЮЧЕНИЕ

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

ЛИТЕРАТУРА

1. Барлоу Р., Прошан Ф. Математическая теория надежности. — М.: Советское радио, 1969. — 488 с.

2. Северцев Н.А. Надежность сложных систем в эксплуатации и отработке. — М.: Высшая школа, 1989. — 432 с.

3. Jardine A.K.S., Tsang A.H.C. Maintenance, Replacement, and Reliability: Theory and Applications. — London, N.-Y.: Boca Raton, CRC/Taylor & Francis, 2006. — 330 p.

4. Blischke W, Murthy D.N.P. Warranty Cost Analysis. — N.-Y.: CRC Press, 1994.

5. Parsa H., Jin M. An improved approximation for the renewal function and its integral with an application in twoechelon // International Journal of Production Economics. — 2013. — Vol. 146, N 1. — P. 142—152.

6. Frees E.W. Warranty analysis and renewal function estimation // Naval Research Logistics Quarterly. — 1986. — Vol. 33. — P. 361—372.

7. Кокс Д.Р., Смит В.Л. Теория восстановлений. — М.: Советское радио, 1967. — 292 с.

8. Andronov А., Vishnevsky V, Rykov V. Reliability of a system with two parallel renewable channels // The 10-th Intern. Conf. on Mathematical Methods in Reliability (MMR-2017) July 3 — July 5, 2017, Grenoble, France. — P. 58.

9. Andronov А., Vishnevsky V, Rykov V. On Reliability Function of a Parallel System with Three Renewable Components // Материалы междунар. науч. конф. «Аналитические и вычислительные методы в теории вероятностей и ее приложениях (АВМТВ-2017)», Москва, РУДН, 23—27 окт. 2017 г. / под общ. ред. А.В. Лебедева. — М., 2017. — С. 430—434.

10. Grigoriev L., Kucheryavy V, Rusev V, Sedyh I. Formation of estimates of reliability indicators for active elements in gas transport systems on the basis of refusals statistics // Journal of Polish Safety and Reliability Association. — 2014. — Vol. 5, N 2. — P. 41—47.

11. ГОСТ Р 27.606—2013. Надежность в технике. Управление надежностью. Техническое обслуживание, ориентированное на безотказность.— М.: Стандартинформ, 2014. — 38 c.

12. Smith W.L., Leadbetter M.R. On the Renewal Function for the Weibull Distribution // Technometrics. — 1963. — N 5. — P. 393—396.

13. Constantine A.G., Robinson N.I. The Weibull renewal function for moderate to large arguments. // Computational Statistics & Data Analysis. — 1997. — N 24. — P. 9—27.

14. Kambo N.S., Rangan A. & Ehsan Moghimi Hadji E. Moments based approximation to the renewal function // Communications in Statistics. — Theory and Methods. — 2012. — Vol. 41. — P. 851—868.

15. Cui L., Xie M. Some Normal Approximations for Renewal Function of Large Weibull Shape Parameter // Communications in Statistics. — 2003. — Vol. 32, N 1. — P. 1—16.

16. Smeitink E., Dekker R. Simple approximation to the renewal function // IEEE Trans. on Reliability. — 1990. — Vol. 39, N 1. — P. 71—75.

17. Maghsoodloo S., Helvaci D. Renewal and Renewal-Intensity Functions with Minimal Repair // Journal of Quality and Reliability Engineering. — 2014. — ID 857437. — 10 p.

18. Русев В.Н., Скориков А.В. Анализ элементов систем газоснабжения с помощью метода производящих функций моментов // Тр. РГУ нефти и газа им. И.М. Губкина. — 2016. — № 1 (282). — C. 68—79.

19. Rusev V, Skorikov A. On solution of renewal equation in the weibull model // Reliability: Theory & Applications. — 2017. — Vol. 12, N 4 (47). — P. 60—67.

20. Deligonul Z.S., Bilgen S. Solution of the Volterra equation of renewal theory with Galerkin technique using cubic splines // Journal of Statistical Computation and Simulation. — 1984. — N 20. — P. 37—45.

21. Xie M. On the solution of renewal-type integral equations // Communications in Statistics. — Simulation and Computation. — 1989. — Vol. 18 (1). — P. 281—293.

22. Xie M., Preuss W., Cui L. Error analysis of some integration procedures for renewal equation and convolution integrals // Journal of Statistical Computation and Simulation. — 2003. — Vol. 73 (1). — P. 59—70.

23. Tortorella M. Numerical solutions of renewal-type integral equation // Informs Journal of Computing. — 2005. — Vol. 17 (1). — P. 66—74.

24. Канторович Л.В., Акилов Г.П. Функциональный анализ: 3-е изд., перераб. — М.: Наука, 1984. — 752 с.

25. Лубков Н.В. Оптимизация срока службы оборудования по критерию эксплуатационных затрат // Автоматика и телемеханика. — 2008. — № 5. — С. 180—190.

26. Hanscom M.A., Cleroux R. The block replacement problem // J. Statist. Comput. Simul. — 1978. — Vol. 3. — P. 233—248.

27. Berg M. A marginal cost analysis for preventive maintenance policies // European Journal of Operational Research. — 1980. — N 4. — P. 136—142.

Статья представлена к публикации членом редколлегии В.М. Вишневским.

Русев Владимир Николаевич — ст. преподаватель, И rusev.v@gubkin.ru,

Скориков Александр Васильевич — канд. физ.-мат. наук, доцент, И skorikov.a@gubkin.ru,

Российский государственный университет нефти и газа (национальный исследовательский университет) имени И.М. Губкина.