Анализ работы восстанавливаемого элемента аналитическими методами и с помощью компьютерного моделирования

Санкт-Петербургский Государственный Университет Информационных Технологий, Механики и Оптики

Отчет по лабораторной работе

«Анализ работы восстанавливаемого элемента аналитическими методами и с помощью компьютерного моделирования»

Выполнили:            

Санкт-Петербург, 2010

Анализ надежности элементов систем по результатам испытаний

Цель работы.

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

Исходные данные.

Модель отказов. Модель процесса восстановления.

Экспериментальные данные.

1.  Генерация случайной величины времени наработки на отказ.

Равномерное распределение на интервале:

Листинг кода в среде MatLab:

%Параметры равномерного распределения на интервале

dWork = 3.3589;

aWork = 8.3875;

a = sqrt(3*dWork)+aWork;

b = 2*aWork-a;

fprintf('Теоретическое математическое ожидание времени безотказной работы: %4.2f\n' , aWork);

fprintf('Теоретическая дисперсия времени безотказной работы: %4.2f\n', dWork);

fprintf('Параметр a равномерного распределения на интервале: %4.2f\n', a);

fprintf('Параметр b равномерного распределения на интервале: %4.2f\n', b);

%Равномерное распределение на интревале

Y=getMyNorm(a, b, 0, 1);

function res = getMyNorm(a, b, be, en);

U = rand(be,en);

U = a+(b-a)*U;

res = U;

%Считаем среднее время наработки на отказ

avgWork=getMyAvg(Y);

fprintf('Эмпирическое мат.ожидание времени безотказной работы: %4.2f\n', avgWork);

%Считаем дисперсию наработки на отказ

dispWork = getMyDisp(Y);

fprintf('Эмпирическая дисперсия времени безотказной работы: %4.2f\n', dispWork);

Получены значения:

2.   Генерация случайной величины времени восстановления.

Распределение Вейбулла:

Листинг кода в среде MatLab:

%Параметры распределения Вейбулла

dWork2 = dWork/2;

aWork2 = aWork/2;

alfa = 4.53;

beta = 4.58;

fprintf('Теоретическое математическое ожидание времени восстановления: %4.2f\n', aWork2);

fprintf('Теоретическая дисперсия времени восстановления: %4.2f\n', dWork2);

fprintf('Параметр alfa распределения Вейбулла: %4.2f\n', alfa);

fprintf('Параметр beta распределения Вейбулла: %4.2f\n', beta);

%Распределение Вейбулла

Z=getMyWbl(alfa, beta, 1, 1000);

function res = getMyWbl(alfa, beta, be, en);

U=rand(be,en);

U=wblinv(U,alfa,beta);

res=U;

%Считаем среднее время восстановления

avgWait=getMyAvg(Z);

fprintf('Эмпирическое математическое ожидание времени восстановления: %4.2f\n', avgWait);

%Считаем дисперсию восстановления

dispWait = getMyDisp(Z);

fprintf('Эмпирическая дисперсия времени восстановления: %4.2f\n', dispWait);

Получены значения:

3. Вывод гистограмм распределений:

Функция распределения времени безотказной работы (эмпирическая):

Функция распределения времени безотказной работы (теоретическая):

Плотность распределения времени безотказной работы:

Функция распределения времени восстановления (эмпирическая):

Функция распределения времени восстановления (теоретическая):

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

4.  Построение аналитической и экспериментальной характеристики готовности системы.

Листинг кода в среде MatLab:

%Теоретический коэффициент готовности

KgTeor = avgWork/(avgWork+avgWait);

fprintf('теоретический коэффициент готовности: %4.2f\n', KgTeor);

%Считаем время

time=floor(sum(Y)+sum(Z));

timeNab=round(4/5*time);

fprintf('Среднее время 1000 итераций устройства: %4.2f\n', time);

fprintf('Время наблюдения: %4.2f\n', timeNab);

%Считаем эмпирический коэффициент готовности

Kg=zeros(1,time);

Y=getMyNorm(a, b, items, 1200);

Z=getMyWbl(alfa, beta, items, 1200);

for t=1:time

work=0;

for item=1:items

if(isReady(Y, Z(item, 1:length(Z)), t))

work=work+1;

end;

end;

Kg(t)=work/items;

end;

figure;

axis([0 time 0 1]);

plot(1:time, Kg, 1:time, KgTeor);

title('Зависимость коэффициента готовности от времени');

figure;

plot(1:100, Kg(1:100), 1:100, KgTeor);

KgT(1)=Kg(timeNab);

for i=1:9

Y=getMyNorm(a, b, items, 1100);

Z=getMyWbl(alfa, beta, items, 1100);

work=0;

for item=1:items

if(isReady(Y, Z(item, 1:length(Z)), timeNab))

work=work+1;

end;

end;

KgT(i+1)=work/items;

end;

Получаем значение:

Kg = 0.67

Графики зависимости коэффициента готовности от времени:

5.  Построение доверительного интервала для оценки готовности системы.

Листинг кода в среде MatLab:

avgKg=getMyAvg(KgT);

dispKg=getMyDisp(KgT);

%Считаем доверительный интервал коэффициента готовности:

delta = 1.645* sqrt(dispKg)/ sqrt(length(KgT));

fprintf('Доверительный интервал с надежностью 0.95 для коэффициента готовности: %4.8f - %4.8f \n', avgKg-delta, avgKg+delta);

Получаем значение:

0.64893901 <  Kg <  0.72106099

6.  Определение вероятности того, что реальное значение Kg < Kg/2.

Листинг кода в среде MatLab:

%Оцениваем вероятность Kg(t) < Kg/2

pKg1=0;

for i=1:length(Kg)

if Kg(i)<KgTeor/2

pKg1 = pKg1 + Kg(i);

end

end;

pKg1 = pKg1/length(Kg);

pKg2 = avgKg/KgTeor;

fprintf('Оценка вероятности P(Kg(t) < Kg/2):\n теоретическая: 0\n эмпирическая: %4.8f\n\n', pKg1);

Получаем значение:

P теоретическая = 0

Р практическая = 0.00022934