Bartłomiej Piekarski 171160 Data utworzenia: 01.06.2010r.
A
ukasz Tkacz 171032
A
ukasz Przywarty 171018
Zadanie projektowe:
Niezawodność i diagnostyka układów cyfrowych
Temat:  Ocena niezawodności systemu pomiarowego typu '2z3' z zawod-
nym układem diagnostyki
ProwadzÄ…cy:
dr inż. Kazimierz Kapłon
Strona 1
Spis treści
I. Szczegółowy opis zadania projektowego . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .3
II. Koncepcja symulacji . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3
a) Wykresy czasowe . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4
b) Tabela stanów . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5
c) Graf stanów . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .6
d) Schemat blokowy . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7
III. Implementacja programowa modelu symulacyjnego w środowisku Matlab . . . . . . . . . . . . 9
a) Listing kodu. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10
b) Wyniki symulacji. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16
IV. Szacowanie charakterystyk niezawodnościowych . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17
V. Podsumowanie i wnioski. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 24
VI. Bibliografia . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 25
Strona 2
I. Szczegółowy opis zadania projektowego
W projekcie należy określić charakterystyki niezawodnościowe systemu pomiarowego,
który składa się z trzech elementów pomiarowych (czujników) i układu diagnostyki i przełączania.
Czujniki są funkcjonalne i niezawodnościowo takie same. Do poprawnej pracy systemu
pomiarowego konieczna jest sprawność co najmniej dwóch czujników. System wykonuje swoją
pracę w sposób ciągły. W ujęciu niezawodnościowym wszystkie elementy są zawodne i napra-
wialne. Znane są: rozkład prawdopodobieństwa czasu poprawnej pracy czujnika do uszkodzenia
F (t) i układu diagnostyki F (t) oraz rozkłady G (t) i odpowiednio G (t) czasu do zakończenia
c d c d
naprawy. Elementy systemu sÄ… naprawiane przez K ekip remontowych. Ponadto:
" układ diagnostyki i przełączania czujników w stanie uszkodzenia nie zakłóca pracy
czujników (tzn. następuje utrata zdolności do wykrywania uszkodzeń i przełączania
elementów, ale jeśli co najmniej 2 z 3 czujników są w stanie sprawności, to system nadal
poprawnie wykonuje pomiary),
" czasy do uszkodzenia elementu Ä , Ä majÄ… rozkÅ‚ady wykÅ‚adnicze o parametrach  , 
c d c d
" czasy naprawy Åš , Åš majÄ… rozkÅ‚ady normalne obciÄ™te (Åš > 0) o parametrach ź , à , ź , Ã
c d c c d d
(wartości dobrać indywidualnie)
W ramach projektu należy:
" opracować model symulacyjny działania systemu uwzględniający występowanie uszkodzeń
i napraw elementów (przedstawić ten model w postaci schematu blokowego)
" opracować implementację programową modelu symulacyjnego w środowisku Matlab,
" dla wybranych wartości parametrów wykonać na modelu niezbędną liczbę eksperymentów
" na podstawie otrzymanych wyników oszacować charakterystyki niezawodności badanego
systemu
II. Koncepcja symulacji
Zadanie polega na utworzeniu modelu symulacyjnego systemu '2z3' z zawodnym układem
diagnostyki. Szczegółowe założenia projektowe zawarte są w punkcie I. niniejszej pracy.
Na samym początku należy przyjrzeć się zachowaniu systemu w poszczególnych chwilach czasu.
Zakładając, że dostępny jest jeden konserwator, w jednej chwili czasu może on naprawiać tylko
jeden element systemu. W momencie gdy uszkodzi się inny element, a konserwator będzie zajęty,
element pozostaje uszkodzony i czeka na naprawę. Szczególnie traktujemy układ diagnostyki:
naprawiamy go w pierwszej kolejności (gdy jest dostępny konserwator, który taką naprawę może
Strona 3
wykonać. W przeciwnym razie układ naprawiany jest zaraz po ukończeniu naprawy innego
elementu).
a) Wykresy czasowe
Przedstawimy teraz przykładowe wykresy czasowe, zakładając, że korzystamy z usług jednego
konserwatora a stany elementów systemu oznaczone są następująco:
" 0  sprawny
" 1 - uszkodzony w naprawie
" 2  uszkodzony bez naprawy
Na rysunku poniżej c c c oznaczają kolejne czujniki, natomiast D  układ diagnostyki
1, 2, 3,
Rysunek 1: Wykresy czasowe dla rozpatrywanego systemu
Jak można zauważyć na Rys. 1 uzyskaliśmy wstępne wyniki symulacji dla konkretnych zdarzeń
Strona 4
(uszkodzeÅ„ i napraw). OznaczyliÅ›my literÄ… Q czasy naprawy, natomiast literÄ… Ä  czasy miÄ™dzy
uszkodzeniami.
b) Tabela stanów
Pomocne będzie rozpatrzenie wszystkich kombinacji stanów poszczególnych elementów wraz
z uwzględnieniem stanu całego systemu. Niestety, wiele z kombinacji jest niemożliwych, tzn., że
nie spełniają warunków określonych w treści zadania projektowego. Kombinacje stanów
przedstawimy w postaci tabeli. Kolorem szarym zaznaczono niemożliwe kombinacje stanów
poszczególnych elementów., wyróżniono natomiast stany możliwe. Zakładamy, że K = 1.
Stan Stan
Stan S Stan elementu Stan S Stan elementu
sys sys
systemu systemu
D S1 S2 S3 D S1 S2 S3
S 0 0 0 0 0 S 0 0 1 0 0
0000 0010
S 0 0 0 1 0 S 0 0 1 1
0001 0011
S 0 0 0 2 S 0 0 1 2 1
0002 0012
S 0 0 2 0 S 1 1 2 1
0020 1121
S 0 0 2 1 1 S 1 1 2 2
0021 1122
S 0 0 2 2 S 1 2 0 0 0
0022 1200
S 0 1 0 0 0 S 1 2 0 1
0100 1201
S 0 1 0 1 S 1 2 0 2 1
0101 1202
S 0 1 0 2 1 S 1 2 1 0
0102 1210
S 0 1 1 0 S 1 2 1 1
0110 1211
S 0 1 1 1 S 1 2 1 2
0111 1212
S 0 1 1 2 S 1 2 2 0 1
0112 1220
S 0 1 2 0 1 S 1 2 2 1
0120 1221
S 0 1 2 1 S 1 2 2 2 1
0121 1222
S 0 1 2 2 1 S 2 0 0 0
0122 2000
S 0 2 0 0 S 2 0 0 1 0
0200 2001
S 0 2 0 1 1 S 2 0 0 2
0201 2002
S 0 2 0 2 S 2 0 1 0 0
0202 2010
S 0 2 1 0 1 S 2 0 1 1
0210 2011
S 0 2 1 1 S 2 0 1 2 1
0211 2012
S 0 2 1 2 1 S 2 0 2 0
0212 2020
S 0 2 2 0 S 2 0 2 1 1
0220 2021
Strona 5
Stan Stan
Stan S Stan elementu Stan S Stan elementu
sys sys
systemu systemu
S 0 2 2 1 1 S 2 0 2 2
0221 2022
S 0 2 2 2 S 2 1 0 0 0
0222 2100
S 1 0 0 0 0 S 2 1 0 1
1000 2101
S 1 0 0 1 S 2 1 0 2 1
1001 2102
S 1 0 0 2 0 S 2 1 1 0
1002 2110
S 1 0 1 0 S 2 1 1 1
1010 2111
S 1 0 1 1 S 2 1 1 2
1011 2112
S 1 0 1 2 S 2 1 2 0 1
1012 2120
S 1 0 2 0 0 S 2 1 2 1
1020 2121
S 1 0 2 1 S 2 1 2 2 1
1021 2122
S 1 0 2 2 1 S 2 2 0 0
1022 2200
S 1 1 0 0 S 2 2 0 1 1
1100 2201
S 1 1 0 1 S 2 2 0 2
1101 2202
S 1 1 0 2 S 2 2 1 0 1
1102 2210
S 1 1 1 0 S 2 2 1 1
1110 2211
S 1 1 1 1 S 2 2 1 2 1
1111 2212
S 1 1 1 2 S 2 2 2 0
1112 2220
S 1 1 2 0 S 2 2 2 1 1
1120 2221
S 2 2 2 2
2222
Tabela 1: Tabela stanów systemu
Okazuje się, że większość kombinacji jest niepoprawna. Możliwe zestawienia stanów elementów są
mniej liczne. Nasz system, przy założeniu, że korzystamy z usług jednego konserwatora przez
większa część czasu pozostaje niesprawny. Stan sprawności systemu zachodzi tylko w 11
przypadkach.
c) Graf stanów
Aby lepiej zobrazować zachowanie systemu w momencie naprawy lub uszkodzenia elementu
przedstawimy graf stanów. Graf stanów przedstawia tylko możliwe stany systemu, oznaczone
kółkiem, w którym widnieje nazwa stanu. Poszczególne przejścia między stanami oznaczone są
strzałkami. Adnotacje przy strzałkach (n lub u ) oznaczają naprawę lub uszkodzenie
i i
poszczególnych elementów. Stan sprawności całego systemu wyróżniono pogrubieniem stanów.
Stany niesprawności zaznaczono przerywaną linią.
Strona 6
Rysunek 2: Graf stanów
Zależności między stanami są dość skomplikowane, ale jednoznaczne. W momencie uszkodzenia
lub naprawy elementu przechodzimy do kolejnego stanu. Procedura ta powtarza się przez cały czas
trwania symulacji.
c) Schemat blokowy algorytmu symulacji
W tym miejscu możemy skupić się na algorytmie symulacji w postaci schematu blokowego.
Schemat blokowy zawiera wyczerpujące dane na temat działania systemu.
Strona 7
Rysunek 3: Model symulacyjny systemu - schemat blokowy
Strona 8
Legenda:
" F (t) oraz F (t) - rozkład prawdopodobieństwa czasu poprawnej pracy do uszkodzenia
c d
czujnika oraz układu diagnostyki
" G (t) oraz G (t) - rozkład prawdopodobieństwa czasu naprawy czujnika i układu
c d
diagnostyki
" K  liczba konserwatorów
" T  całkowity czas symulacji
sym
" L  liczba symulacji
sym
" T , T , T oraz T  moment uszkodzenia poszczególnych czujników oraz układu
1 2 3 d
diagnostyki
" T  czas najbliższego uszkodzenia (wydarzenia)
x
" S oraz S  stany czujników i układu diagnostyki
x d
" L
0  liczba oczekujÄ…cych na naprawÄ™
Opiszemy pokrótce algorytm działania systemu. Na samym początku określamy niezbędne
parametry, takie jak rozkłady prawdopodobieństw, liczba symulacji itd. Losujemy czasy poprawnej
pracy i wybieramy minimum z tych czasów. W momencie gdy uszkodzeniu uległ konkretny
element systemu naprawiamy go, zakładając, ze są dostępni konserwatorzy. Priorytet (jeśli chodzi
o kolejność wykonywanych napraw) ma układ diagnostyki. Czas naprawy jest również parametrem
losowym. Gdy nie ma dostępnych konserwatorów zwiększamy liczbę oczekujących na naprawę
i staramy się o rozpoczęcie naprawy w innym momencie czasu. Jeśli wszystkie elementy zostały
naprawione zwiększamy czas od startu symulacji i zmniejszamy czasy uszkodzeń. Wyniki
poszczególnych prób zapisujemy w postaci plików wyjściowych o konkretnym numerze.
Wszystkie kroki powtarzamy do momentu aż przekroczymy maksymalny czas symulacji
i maksymalną ilość symulacji.
Dokładny opis działania symulatora znajduje się w punkcie IV niniejszej pracy.
IV. Implementacja programowa modelu symulacyjnego w środowisku Matlab
Program symulacyjny oprócz tego, że generuje losowe zdarzenia (uszkodzenia, naprawy)
będzie także obliczał następujące miary niezawodnościowe:
" TTF  czas do uszkodzenia,
" TTR  czas do naprawy,
" MTBF  średni czas między uszkodzeniami,
Strona 9
" MTTF  średni czas do uszkodzenia,
" MTTR  średni czas do naprawy.
a) Listing kodu
Opiszemy działanie symulatora na podstawie dołączonych fragmentów kodu. Na samym początku
ustalamy maksymalny czas symulacji oraz liczbę tych symulacji. Dobieramy również rozkłady
prawdopodobieństwa dla czasów poprawnej pracy i czasów naprawy.
clear
Tsym = 500;
Lsym = 1000;
F = 13;
G = 10;
Listing 1: Podstawowe parametry
Ustalamy ścieżki zapisu plików zawierających czasy do uszkodzenia  TTF oraz czasy do naprawy
TTR. W tych plikach przechowywać będziemy obliczone w dalszej części programu wartości.
sciezka_ttf = sprintf('wyniki\\ttf.txt');
[file_ttf, message] = fopen(sciezka_ttf,'w');
%Sprawdzenie poprawności otwarcia
if file_ttf == -1
disp(wiadomosc)
return;
end;
sciezka_ttr = sprintf('wyniki\\ttr.txt');
[file_ttr, wiadomosc] = fopen(sciezka_ttr,'w');
if file_ttr == -1
disp(wiadomosc)
return;
end;
Listing 2: Ścieżki zapisu plików ttf i ttr
W tym miejscu możemy zaimplementować schemat blokowy działania systemu  w postaci kilku
zagnieżdżonych pętli 'for'. Definiujemy zmienne i losujemy momenty uszkodzenia
for L = 1:Lsym
%Inicjalizacja zmiennych.
K = 1; %liczba konserwatorów
Strona 10
 T = 0; %czas symulacji
tx = 0; %bieżąca chwila (minimum z wylosowanych momentów)
LO = 0; %liczba oczekujÄ…cych na naprawÄ™
S = [0,0,0,0]; %stany poszczególnych elementów systemu - pierwszy układ
%diagnostyki, następnie czujniki
Tx = [1,2,3,10]; %czasy uszkodzenia (najbliższego wydarzenia)
%Losowanie momentów uszkodzenia
for i = 1:4
Tx(i) = exprnd(F);
end
Listing 3: Zmienne i momenty uszkodzenia
Sprawdzimy również czy system jest sprawny. Wykorzystamy do tego celu dodatkową zmienną
sprawne, która będzie przechowywała liczbę elementów sprawnych.
%Zmienna przechowująca ilość elementów sprawnych
sprawne = 0;
for i = 2:4
%Jeśli stan elementu jest równy 0 (jest sprawny) zwiększamy
%zmiennÄ… sprawne o 1
if S(i) == 0
sprawne = sprawne + 1;
end
end
%W przypadku gdy liczba czujników sprawnych jest równa przynajmniej 2
%cały system jest sprawny (zasada '2z3')
if sprawne > 1
stan = 0;
%w przeciwnym wypadku układ jest niesprawny
else
stan = 1;
end
Listing 4: Sprawdzenie stanu systemu
Mając stany poszczególnych elementów jak i całego systemu możemy zapisać wyniki do pliku.
Definiujemy ścieżkę zapisu pliku i wpisujemy do niego wartości:
Strona 11
 %Ustalenie ścieżki zapisu dla plików wynikowych symulacji
sciezka_proba = sprintf('wyniki\\sym%d.txt',L);
[file_proba, wiadomosc] = fopen(sciezka_proba,'w');
if file_proba == -1
disp(wiadomosc)
return;
end;
%Zapisujemy wyniki do pliku w postaci: stan systemu, stan układu
%diagnostyki, stany czujników, czas od momentu startu symulacji.
fprintf(file_proba,'%d\t%d\t%d\t%d\t%d\t%f\r\n', stan, S(1), S(2), S(3),
S(4),T);
Listing 5: Zapis stanów elementów i systemu do pliku
Dochodzimy do najważniejszej części programu  symulacji. Dopóki czas jest mniejszy od
maksymalnego czasu symulacji wykonujemy procedury zobrazowane w schemacie blokowym.
Rozpoczynamy od wyszukania minimum z czasów uszkodzeń.
while T < Tsym
%zwiększamy bieżącą chwilę
tx = 2 * Tsym;
%Szukamy minimum z czasów uszkodzeń
for i = 1:4
if Tx(i) < tx
x = i;
tx = Tx(i);
end
end
Listing 6: Szukanie minimum z czasów uszkodzeń
Następnie sprawdzamy czy urządzenie w chwili tx_ było sprawne. Jeśli tak  sprawdzamy czy
układ diagnostyki był sprawny i czy liczba konserwatorów jest większa od 0. W przypadku
uzyskania twierdzącej odpowiedzi wiemy, że urządzenie uległo uszkodzeniu, musimy je więc
naprawić. Zmniejszamy liczbę dostępnych konserwatorów i losujemy czas naprawy. Zaznaczamy
stan elementu jako 1 (w naprawie). W przeciwnym razie zwiększamy liczbę oczekujących na
naprawÄ™ i oznaczamy stan elementu jako 2 (czekajÄ…cy na naprawÄ™).
if S(x) == 0
if S(1) == 0 && K>0
Strona 12
K = K-1;
Tx(x) = exprnd(G) + tx;
S(x) = 1;
else
Tx(x) = 5 * Tsym;
LO = LO + 1;
S(x) = 2;
end
Listing 7: Sprawdzenie sprawności w chwili tx_
Jeśli element w chwili tx_ był uszkodzony losujemy czas poprawnej pracy, zwiększamy liczbę
konserwatorów i oznaczamy stan elementu jako 0 (sprawny)
else
K = K+1;
Tx(x) = exprnd(F) + tx;
S(x) = 0;
end
Listing 8: Sprawdzenie sprawności w chwili tx_ c.d.
Szczególnie traktujemy układ diagnostyki. Jeśli układ ten oczekuje na naprawę i są dostępni
konserwatorzy naprawiamy go w pierwszej kolejności.
if S(1) == 2 && K > 0
K = K - 1;
Tx(1) = exprnd(G) + tx;
S(1) = 1;
LO = LO - 1;
end
Listing 8: Próba naprawy układu diagnostyki
Nie tylko układ diagnostyki ulega uszkodzeniu. Podobnie dzieje się z czujnikami, dopóki liczba
oczekujących na naprawę czujników jest większa od 0, układ diagnostyki jest sprawny i są dostępni
konserwatorzy, naprawiamy poszczególne elementy w kolejności 1-wszy, 2-gi, 3-ci czujnik.
i = 2;
while LO > 0 && S(1) == 0 && K > 0
if S(i) == 2
K = K - 1;
Tx(i) = exprnd(G) + tx;
S(i) = 1;
Strona 13
LO = LO - 1;
end
i = i + 1;
end
Listing 9: Naprawa poszczególnych czujników systemu
Możemy teraz zwiększyć czas od startu systemu a także zmniejszyć czasy najbliższego wydarzenia
i kolejny raz policzyć ile elementów jest sprawnych. Określenie liczby sprawnych elementów
systemu wygląda identycznie jak w Listingu 4, więc pominiemy ten fragment kodu. Wyniki
zapisujemy do pliku w sposób zaprezentowany w końcowej części Listingu 5.
%Zwiększamy czas od startu systemu
T = T + tx;
%Zmniejszamy czasy najbliższego wydarzenia
for i = 1:4
Tx(i) = Tx(i) - tx;
end
Listing 10: Korekta czasów
Nasz symulator wygenerował do tej pory uszkodzenia i naprawy a także czasy tych zdarzeń. Mając
do dyspozycji te dane jesteśmy w stanie obliczyć parametry TTF oraz TTR. Rezultaty zapisujemy
do plików tekstowych.
%Otwieramy plik z wynikami konkretnej próby
[file_proba, wiadomosc] = fopen(sciezka_proba,'r');
if file_proba == -1
disp(wiadomosc)
return;
end
%Z pliku zczytujemy wartości do tablicy
tab = fscanf(file_proba,'%d%d%d%d%d%f',[ 6 inf ]);
fclose(file_proba);
tab = tab';
rozmiar = size(tab,1);
ostatniStan = tab(1,1); %Pierwszy stan systemu
ostatniRazSprawny = 0; %Moment ostatniej sprawności
ostatniRazUszkodzony = 0;%Moment ostatniego uszkodzenia
%W tym miejscu obliczymy poszczególne czasy do uszkodzenia i czasy do
Strona 14
 %naprawy - TTF oraz TTR
for i = 2:rozmiar
%Zczytujemy czas od momentu startu symulacji
czas = tab(i,6);
%oraz stan systemu
stan = tab(i,1);
%Jeśli ostatni odczytany stan różni się od bieżącego sprawdzamy czy
%jest równy 1 jeśli tak zapisujemy do pliku obliczoną wartość TTF.
%W przeciwnym razie - wartość TTR
if ostatniStan ~= stan
if stan == 1
czasSprawnego = (czas - ostatniRazSprawny);
ostatniRazUszkodzony = czas;
fprintf(file_ttf,'%f\r\n', czasSprawnego);
else
czasUszkodzonego = (czas - ostatniRazUszkodzony);
ostatniRazSprawny = czas;
fprintf(file_ttr,'%f\r\n', czasUszkodzonego);
end
end
ostatniStan = stan;
end
end
fclose(file_ttf);
fclose(file_ttr);
Listing 11: Obliczenie parametrów TTF i TTR
Aby umożliwić póz
niejszą analizę wyników umieścimy wartości TTF oraz TTR w tablicach:
%Otwieramy pliki z wartościami TTF oraz TTR
sciezka_ttf = sprintf('wyniki\\ttf.txt');
[file_ttf, wiadomosc] = fopen(sciezka_ttf,'r');
if file_ttf == -1
disp(wiadomosc)
return;
end;
sciezka_ttr = sprintf('wyniki\\ttr.txt');
[file_ttr, wiadomosc] = fopen(sciezka_ttr,'r');
if file_ttr == -1
disp(wiadomosc)
Strona 15
 return;
end;
%Zapisujemy wartości TTF oraz TTR w tablicy
ttf = fscanf(file_ttf,'%f', [ 1 inf ]);
ttf = ttf';
ttr = fscanf(file_ttr,'%f',[ 1 inf ]);
ttr = ttr';
fclose(file_ttf);
fclose(file_ttr);
Listing 12: Zapis wartości TTF i TTR w tablicach
Na samym końcu obliczamy pozostałe parametry MTBF, MTTF, MTTR
fprintf('MTBF: %f\n',mean(ttf) + mean(ttr));
fprintf('MTTF: %f\n',mean(ttf));
fprintf('MTTR: %f\n',mean(ttr));
Listing 13: Obliczenie pozostałych parametrów niezawodnościowych.
b) Wyniki symulacji
Uruchamiając program uzyskaliśmy pliki wynikowe, których fragmenty przedstawiono w Tabeli 2.
Stan Stan Stan 1. Stan 2. Stan 3. Czas od TBF MTBF MTTF MTTR
układu diagno- czujni- czujni- czujni- początku
styki ka ka ka symulacji
0 0 0 0 0 0.000000 -
0 1 0 0 0 7.687255 -
0 0 0 0 0 8.781624 -
0 0 0 1 0 9.834445 -
0 2 0 1 0 9.864509 -
1 2 2 1 0 17.490966 17.490966
0 1 2 0 0 29.341057 -
0 0 1 0 0 32.566749 -
36.870919 8.776564 28.094354
1 0 1 0 2 37.286037 19.795071
1 0 1 2 2 37.961525 0.675488
1 0 0 1 2 44.649317 6.687792
0 0 0 0 1 55.324184 -
1 0 2 0 1 56.910808 12.261491
1 2 2 0 1 58.083068 1.172260
1 2 2 2 1 66.521111 8.438043
1 1 2 2 0 70.361735 3.840624
Tabela 2: Fragment plików wynikowych
Strona 16
IV. Szacowanie charakterystyk niezawodnościowych
Podsumujmy nasze dotychczasowe działania: analizując treść zadania projektowego
utworzyliśmy schemat blokowy działania systemu. Na tej podstawie napisaliśmy program
w środowisku Matlab. Program symuluje działanie systemu, zapisuje wyniki do poszczególnych
plików i oblicza wartości TTF oraz TTR. Wykorzystamy je do szacowania charakterystyk
niezawodnościowych.
Aby ukazać z jaką częstością pojawiają się kolejne wartości TTF i TTR wykorzystamy histogramy.
Histogram w programie Matlab można utworzyć wydając polecenie: hist(dane). Jeśli
w miejsce danych podstawimy ttf uzyskamy wykres:
Rysunek 4: Histogram częstości wystąpień wartości TTF
Postępując podobnie w przypadku ttr wyświetlimy drugi histogram:
Rysunek 5: Histogram częstości wystąpień wartości TTR
Strona 17
Z Rys. 4 oraz Rys. 5 wynika, że zdecydowanie częściej występują małe wartości TTF oraz TTR. Im
większa wartość tym rzadziej ona występuje. Wiedząc jak wyglądają histogramy wartości TTF
i TTR możemy przystąpić do właściwej części zadania  dopasowywania rozkładów.
W trakcie poszukiwań informacji na temat dopasowywania rozkładów do zmiennych losowych
natrafiliśmy na bardzo pomocne narzędzie zaimplementowane w programie Matlab  dfittool.
Na podstawie danych (w naszym przypadku tablic z wartościami TTF oraz TTR) generuje
histogramy, do których można dopasować funkcje rozkładu np. Weibulla, gamma, normalnego.
Narzędzie korzysta z odrębnego graficznego interfejsu użytkownika, wszystkie operacje
wykonujemy za pomocą dostępnych przycisków. Dfittool automatycznie generuje wykresy
oraz potrzebne parametry np. dla rozkładu Weibulla  parametr skali oraz kształtu (odpowiednio
a i b w programie Matlab). Użycie narzędzia jest bardzo proste wystarczy wpisać (np. w oknie
poleceń): dfittool(dane) w naszym przypadku: dfittool(ttf) lub dfittool(ttr)
Po uruchomieniu narzędzia otwiera się osobne okno przedstawione na Rys. 6. Okno to zawiera
kilka ważnych elementów.
Najistotniejsze dla nas sÄ…
przyciski znajdujÄ…ce siÄ™ nad
wykresem. Aby zaimportować
zbiór danych należy skorzystać
z przycisku Data. Kolejny
raz wyświetli się nowe okno,
które pozwala wybrać, do któr-
ych zmiennych dopasowywać
będziemy rozkłady. Operacje
dodawania danych akceptu-
jemy przyciskiem Create
data set.
Przycisk New Fit... po-
Rysunek 6: Główne okno narzędzia dfittool
zwala dodać nowe dopaso-
wanie. Podobnie jak w przypadku tworzenia zbioru zmiennych korzystamy z odrębnego okna
wyboru, które umożliwia nazwanie konkretnego dopasowania (Fit name) oraz zdefiniowanie po-
żądanego rozkładu (Distribution). Warto zaznaczyć, że do jednego wykresu możemy dodać
wiele dopasowań, pozwala to łatwiej określić, który rozkład jest najbardziej odpowiedni.
Przyciski Manage fits... służą do zarządzania dopasowaniami, natomiast Evaluate...
Strona 18
oraz Exclude... kolejno do badania konkretnych zmiennych oraz wykluczania wartości,
których nie chcemy obejmować dopasowywaniem.
Ważnym elementem okna narzędzia dfittol jest również pole Display type. Pozwala ono
wybrać sposób wyświetlania wykresów np.Density oznacza wyświetlanie w postaci funkcji
gęstości prawdopodobieństwa natomiast Cumulative probability  funkcji rozkładu
prawdopodobieństwa.
Spróbujemy teraz dopasować histogram wartości TTF oraz TTR do kilku rozkładów. Uznaliśmy,
że najodpowiedniejsze będą rozkłady:
" normalny,
" eksponencjalny,
" Weibulla.
Uzyskaliśmy następujące rezultaty (w postaci funkcji gęstości prawdopodobieństwa):
Rysunek 7: Dopasowania rozkładów do wartości TTF
Z analizy przeprowadzonego doświadczenia (Rys. 7 oraz Rys. 8) wynika , że najlepiej dopasowany
do histogramu jest rozkład eksponencjalny. Rozkład normalny kompletnie nie pasuje, natomiast
rozkład Weibulla osiąga zbyt duże gęstości.
Strona 19
Rysunek 8: Dopasowania rozkładów do wartości TTR
Narzędzie dfittol wygenerowało dla rozkładów następujące parametry:
" dla rozkładu normalnego:
 mu (średnia ź)  8.77656
 sigma (odchylenie standardowe Ã)  9.48085
" dla rozkładu Weibulla
 a (parametr skali ) 8.45563
 b (parametr kształtu k) - 0.92366
" dla rozkładu eksponencjalnego
 mu (parametr )  8.77656
Niestety  tylko na podstawie dopasowań do histogramów nie możemy jednoznacznie określić,
który rozkład najlepiej pasuje do naszych wartości. Przedstawimy więc te same dane i rozkłady na
wykresie funkcji prawdopodobieństwa. W tym celu w polu Display type narzędzia
dfittol wybieramy opcjÄ™ Probability plot.
Strona 20
Rysunek 9: Wykres probabilistyczny dla wartości TTF i rozkładu eksponencjalnego
Rysunek 10: Wykres probabilistyczny dla wartości TTR i rozkładu eksponencjalnego
Widzimy, że o ile dla wartości TTF rozkład eksponencjalny jest dość dobrze dopasowany to dla
wartości TTR  nie (linia czerwona oznacza dopasowanie idealne). Porównamy więc wyniki
z pozostałymi rozkładami.
Strona 21
Rozkład normalny:
Rysunek 11: Wykres probabilistyczny dla wartości TTF i rozkładu normalnego
Rysunek 12: Wykres probabilistyczny dla wartości TTR i rozkładu normalnego
Podobnie jak w przypadku dopasowania do histogramu rozkład normalny odrzucamy, ze względu
na to, że zdecydowanie odbiega od ideału (Rys. 11 oraz Rys 12)
Strona 22
Rozkład Weibulla:
Rysunek 13: Wykres probabilistyczny dla wartości TTF i rozkładu Weibulla
Rysunek 14: Wykres probabilistyczny dla wartości TTR i rozkładu Weibulla
W tym momencie pojawia się problem. Z analizy dopasowań do histogramu wynika, że nasz
rozkład najlepiej przybliża rozkład eksponencjalny. Zaprzeczeniem tego są wykresy
Strona 23
probabilistyczne, na których jednoznacznie widać, że najbardziej dopasowany jest rozkład Weibulla
(Rys. 12 oraz Rys. 13). Uznajemy więc, że rozkład Weibulla będzie najlepszym wyborem.
V. Podsumowanie i wnioski
W ramach projektu opracowaliśmy model symulacyjny systemu uwzględniający
występowanie uszkodzeń i napraw elementów. Model ten został przedstawiony w postaci
schematu blokowego. Na tej podstawie utworzyliśmy implementację programową modelu
w środowisku Matlab. Poszczególne zdarzenia  uszkodzenia i naprawy wygenerowaliśmy
funkcjami programu Matlab. Symulator tworzy pliki wynikowe z rezultatami, które są
następnie opracowywane w postaci histogramów i wykresów a także miar
niezawodnościowych  TTF, TTR, MTBF, MTTF, MTTR. Szacowanie charakterystyk
niezawodnościowych pozwoliło określić, że najlepiej dopasowanym rozkładem jest rozkład
Weibulla. Funkcje rozkładów prawdopodobieństwa liczone są według następujących
wzorów:
" dla rozkładu eksponencjalnego:
f śąxźą=ÁÄ…e-ÁÄ… x , dla x ‡Ä…0
gdzie:  - parametr
f śąxźą=0, dla x"ą0
" dla rozkładu normalnego
-śąx-Âąźą2
1
2 ÈÄ…2
f śą xźą= e , dla x‡Ä…0
gdzie ź  Å›rednia, à  odchylenie standardowe
ÈÄ… ćą 2Ćą
f śąxźą=0, dla x"ą0
" dla rozkładu Weibulla
p-1
f śąxźą=ÁÄ…k śąÁÄ… xźą e-śąÁÄ… xźąk , dla x‡Ä…0
gdzie   parametr skali, k  parametr kształtu
f śąxźą=0, dla x"ą0
Uwaga 1: W pracy często pojawia się pojęcie 'uszkodzenia' systemu. Przypominamy:
system jest niesprawny gdy nie działa więcej niż jeden czujnik lub układ diagnostyki jest
niesprawny. Niesprawność systemu oznaczaliśmy zawsze cyfrą 1, sprawność  0.
Uwaga 2: Kod w postaci pliku Matlab`owego oraz sprawozdanie w wersji elektronicznej sÄ…
zawarte na dołączonej do projektu płycie CD.
Strona 24
VI. Bibliografia
Podczas tworzenia niniejszego projektu korzystaliśmy z następujących z
ródeł:
1. Jarnicki J.  Wykład o estymacji rozkładów zmiennych losowych
2. Praca zbiorowa,  Niezawodność i eksploatacja systemów , pod redakcją W. Zamojskie-
go, Wrocław 1981
3. Pratap R.,  Matlab 7 dla naukowców i inżynierów , Warszawa 2006
Miejsce na płytę CD
Strona 25