Streszczenie
Teoria szarego systemu jako metoda prognozowania stanu obiektu znajduje wiele praktycznych zastosowaĔ, nie tylko w zakresie nauk technicznych, lecz równieĪ w wielu innych gałĊziach nauki, np. naukach społecznych,, przyrodniczych.
Niniejsza publikakcja ma na celu przybliĪenie genezy teorii szarych systemów, ich matematycznej interpretacji. Ponadto wykazano, Īe moteda szarego systemu pierwszego rzĊdu stanowi potĊĪne narzĊdzie w rĊkach inĪyniera, ze wzglĊdu na szerokie zastosowanie tej metody, głównie w zakresie nauk technicznych, czemu poĞwiĊcone jest w głównej mierze niniejsze opracowanie.
Słowa kluczowe: szare systemy, diagnostyka drganiowa, prognozowanie stanu 1. Wprowadzenie
Teoria szarego systemu pierwszego rzĊdu została przedstawiona po raz pierwszy w roku 1982 przez chiĔskiego naukowca – Julonga Denga do modelowania procesów, z których dane mogą byü niepewne i niekompletne, a zasady rządzące zjawiskami nie są do koĔca znane [1] Powstała teoria na długi czas pozostała nieznana poza Chinami, bowiem pierwsza publikacja w jĊzyku angielskim ukazała siĊ dopiero w 1989 r., dziĊki czemu Ğwiat anglojĊzyczny mógł siĊ z nią bliĪej zapoznaü [1]
Metoda ta jako narzĊdzie słuĪące do prognozowania stanu, w odróĪnieniu od innych metod (np. naiwnej, Ğredniej ruchomej) nie zakłada, Īe szereg bĊdzie siĊ zmieniał z okreĞlonym trendem oraz, Īe charakteryzuje siĊ stałą tendencją zmian [3] Istotne natomiast jest, Īe kolejne wyrazy szeregu muszą byü monotonicznie rosnące oraz nieujemne, a iloĞü wprowadzonych danych wejĞciowych musi wynosiü minimum cztery [1]
Teoria szarego systemu moĪe zostaü podzielona na kilka istotnych grup, tj.:
a) Grey Realational Analisys (GRA) – metoda ta słuĪy do okreĞlenia siły związku pomiĊdzy badanymi zmiennymi,[6]
b) Grey Decision Making (GDA) – metoda słuĪąca wspomaganiu podejmowania decyzji,[6] c) Grey Model (GM (1,1) – metoda prognozowania, na podstawie danych wejĞciowych, przy wykorzystaniu równania róĪniczkowego pierwszego rzĊdu z jednym typem wymuszenia. Szczególnym przypadkiem jest metoda GM ze Ğlizgającym siĊ okienkiem (ang. rolling window), która polega na dokonywaniu prognozy dla wąskich przedziałów czasowych, duĪa szerokoĞü przedziału moĪe wpływaü niekorzystnie na dokładnoĞü prognozy. Ten aspekt metody szarego systemu zostanie bliĪej omówiony wraz z moĪliwymi zastosowaniami [6].
2. Matematyczny opis budowy algorytmu szarego systemu pierwszego rzĊdu [2]
Metoda Deng’a, jak kaĪdy algorytm, posiada swoją definicjĊ matematyczną, składającą siĊ zpewnej sekwencji etapów, przedstawionych w kolejnej czĊĞci opracowania. Do przygotowania prognozy uĪyto Ğrodowiska numerycznego Matlab. Na tĊ okolicznoĞü zdefiniowano zapis matematyczny szarego modelu pierwszego rzĊdu ze Ğlizgającym siĊ okienkiem na zapis numeryczny, tak, aĪeby był on właĞciwie interpretowany przez Matlaba. Przed omówieniem zastosowaĔ metody szarego systemu warto opisaü model oraz sposób jego budowy [2]
Szary model, który charakteryzuje zachowanie siĊ systemu na podstawie symptomu, którego postaü okreĞlamy z obserwacji (x(0)(t)), gdzie t to kolejny odczytywany symptom, czyli
t ࣅ <1,2,3,…,),
moĪe byü opisany zwyczajnym równaniem róĪniczkowym k - tego rzĊdu i mieü wymuszenie eͲtego rzĊdu, co umownie moĪna przedstawiü jako GM(k, e), tzn. po prawej stronie stronie równania widzimy k -wymuszeĔ zachowania systemu, charakteryzowane równaniem róĪniczkowym, e – tego rzĊdu, co przyjmuje postaü matematyczną:
(1) 1 (1) 1 0 1 n i e k i n j j i j
d
x
a
b y
dt
− − − = ==
¦
¦
Gdzie: x(1)(k) = ( )0 1( )
t ix
i
=¦
jest zmienną stanu badanego obiektuyj – niezaleĪne wymuszenie, dla poprawnego odwzorowania zachowania badanego obiektu ai , bi – współczynniki wielomianów, wyznaczane z szeregu czasowego x(0)(t), t = 1,2,3…,
Generalnie rzecz ujmując, zwykle przyjmuje siĊ jako punkt wyjĞcia równanie GM(1,1), tzn. równanie pierwszego rzĊdu z jednym tylko rodzajem wymuszenia. Tok postĊpowania w takiej sytuacji moĪe zostaü opisany w postaci kilku etapów, tzn.:
Etap I
OkreĞlenie wektora obserwacji:
x(0) = [x(0) (1), x(0) (2),x(0) (3),…,x(0) (n)],
w przypadku tego wektora warunkiem koniecznym jest, by n 4, gdzie n – liczba obserwacji [1]
Etap II
Utworzenie wektora sum cząstkowych (Accumulating Generating Operation) – AGO: x(1)(t) = ( )0 1
( ),
t ix
i
=¦
t = 1,2,3,…nWidaü, Īe wektor sformułowany w takiej postaci jest monotoniczny oraz rosnący x(1) = [x(1) (1), x(1) (2),x(1) (3),…,x(1) (n)], gdzie warunkiem jest:
x(1) (1) = x(0) (1), [1] Etap III
Dla wyĪej zdefiniowanego wektora AGO dokonujemy opisu szarego modelu pierwszego rzĊdu wg pierwszego schematu GM(1,1) (1) (1)
( )
( )
dx
t
ax
t
u
dt
+
=
Gdzie: a – eksponent wzrostuu – zmienna kontrolna szarego modelu
t – zmienna niezaleĪna – wielkoĞü eksploatacyjna, np. zuĪycie, czas [1] Etap IV
Dokonujemy rozwiązania powyĪszego równania róĪniczkowego z jedynkowym krokiem zmiennej t, w wyniku czego otrzymujemy:
(1)(k+1) = [x(0)(1) – u/a] exp(-ak) + u/a,
jako (1) okreĞliliĞmy prognozĊ naszego rozwiązania wektora sum cząstkowych AGO [1] Etap V
Zamieniamy z równania z etapu III róĪniczkĊ przyrostem skoĔczonym dla t = 1 zapisujemy równania poprzedzające i progresywne, w wyniku zabiegu otrzymując:
równanie poprzedzające: x(1)(k+1) - x(1)(k) + ax(1) (k) = u równanie progresywne: x(1)(k+1) - x(1)(k) + ax(1) (k+1) = u Dokonując kombinacji powyĪszych równaĔ otrzymujemy: x(1)(k+1) - x(1)(k) = -a/2[x(1)(k) + x(1)(k+1)] + u
k = 1,2,3,…,n [1] Etap VI
Rozpisujemy powyĪsze równanie dla kolejnych wartoĞci k i opierając siĊ na wczeĞniejszym wektorze obserwacji x(0) wyznaczamy współczynniki szarego modelu a i u, które dotychczas pozostają nieznane. Wyznaczenie odbywa siĊ wg metody najmniejszych kwadratów, uzyskując macierzowe rozwiązanie, tzn.
[a,u]T = (B TB)-1 BTY
B = ( ) ( ) ( ) 1 (1) 1 (1) 1 (1) [x (1) (2)] 1 [x (2) (3)] 1 ... ... [x ( 1) ( )] 1 x x n x n − + − + − − + Etap VII
Dokonujemy przekształcenia odwrotnego wektora AGO (Inverse AGO) IAGO, uzyskując prognozowaną obserwacje z wektora prognoz sum cząstkowych (AGO), otrzymując:
(0) (k+1) = (1) (k+1) – (1) (k),
przy uĪyciu równania powstałego z kombinacji równania poprzedzającego i progresywnego (etap V) otrzymujemy koĔcową wartoĞü prognozy wg modelu GM(1,1), tzn.:
(0) (k+1) = [x(0) (1) – u/a](e-ak – e-a(k-1) ) k = 2,3,4,…,n [1].
W siódmym kroku algorytmu otrzymuje siĊ cały sens prognozy, przy wykorzystaniu metody Denga [1] Krok drugi wskazuje na własnoĞci wprowadzanych danych – ich monotonicznoĞü, atakĪe by wartoĞci były coraz wiĊksze [2]. Dlatego teĪ uzasadnione jest badanie zjawisk związanych np. ze zuĪywaniem elementów, bowiem zuĪycie elementów maszyn narasta w czasie, niezaleĪnie od tego czy są to elementy proste typu tuleja, tarcza, wałek, czy złoĪone, jak np. łoĪysko. Prognoza stanu obiektu (maszyny) moĪe byü przeprowadzana na podstawie np. obserwowanych symptomów drganiowych [2].
Opisany wyĪej algorytm jest szczególnym przypadkiem metody szarego systemu pierwszego rzĊdu z tzw. Ğlizgającym siĊ okienkiem (ang. rolling window). Jego zastosowanie uwarunkowane jest faktem, Īe duĪa iloĞü danych (obserwacji) wpływa niekorzystnie na błąd prognozy, powodując, Īe staje siĊ ona niedokładna. Dlatego teĪ stosuje siĊ tzw. Ğlizgające siĊ okienko o mniejszej szerokoĞci, które przemieszcza siĊ po wiĊkszej liczbie danych. MetodĊ tą stosuje siĊ kiedy mamy do czynienia z duĪą liczbą danych (zwykle kilkadziesiąt). SłuĪy głównie do przygotowania prognoz krótkookresowych [2].
3 Zastosowanie metody szarego systemu pierwszego rzĊdu do prognozowania jednowymiarowego
3.1. Zastosowanie metody szarego systemu pierwszego rzĊdu do prognozowania zuĪycia łoĪy-ska turbiny samolotu [2]
Jak juĪ wspomniano, opisany i scharakteryzowany powyĪej algorytm ma swoje praktyczne zastosowanie. W tej czesci niniejszej publikacji zwrócono uwage na wykorzystanie metody szarego systemu w diagnozowaniu stanu maszyn. Jako przykład moĪe posłuĪyü łoĪysko turbiny samolotu.
W wyniku zuĪycia dochodzi do jego drgaĔ. Przeprowadzono okreĞloną iloĞü obserwacji parametrow tych drgaĔ, a ĞciĞlej mowiąc prĊdkoĞci drgaĔ, wyraĪonych w (mm/s) w jednakowych odstĊpach – co 20 godzin pracy turbiny. Wyniki przedstawiono w tabeli:
Tabela 1. Zestawienie zaobserwowanych wartoĞci zuĪycia łoĪysk turbiny silnika samolotu (warto-Ğci okreĞlają przyspieszenie drgaĔ w funkcji przebiegu turbiny)
Nr obserwacji 1 2 3 4 5 6 7 8 9 10
WartoĞü (mm/s) 1 3 4 6 8 11 12 14 15 17
Jak widaü wyniki obserwacji mają charakter monotoniczny oraz rosnący, wiec moĪemy tu zduĪym powodzeniem zastosowaü teorie szarego systemu. Dokonujemy zaimplementowania algorytmu w MATLAB i wczytaniu otrzymanych danych. W wyniku tego zabiegu uzyskano nastĊpującą prognozĊ: 0 2 4 6 8 10 12 14 16 18 0 10 20 30 40 50 60 etap obserwacji/prognozy w ar tos c obs er w a c ji/ p rogn oz y [ m m /s ] Pomiary Prognoza
Rysunek 1. Fragment ewolucji łoĪyska turbiny samolotu z 10 obserwacji oraz prognoza jego stanu przy wykorzystaniu metody „rolling window” na kolejnych 8 przegladow – 160 h pracy Jak widaü na powyĪej przedstawionym wykresie, charakter prognozy przypomina nieco wykres funkcji kwadratowej. I tak dzieje siĊ w rzeczywistoĞci, tzn. charakter zuzycia danej czeĞci maszyny przebiega zwykle w sposób nierównomierny, tzn., im wyĪszy stopieĔ zuĪycia, tym jego tempo wzrasta [2]
3.2. Zastosowanie metody szarego systemu pierwszego rzĊdu do prognozowania zuĪycia amortyzatorów w samochodzie osobowym [4]
Zastosowanie metody szarego modelu pierwszego rzĊdu jest w tym przypadku zasadne, po-niewaĪ zuĪycie takiego elementu jak amortyzator jest ze swej natury monotoniczne oraz rosnące wraz z jego eksploatacją. Wyniki uzyskane na podstawie obserwacji przedstawiono w poniĪszej tabeli:
Tabela 2. Zestawienie zaobserwowanych wartoĞci zuĪycia amortyzatorów w samochodzie osobowym w funkcji przebiegu pojazdu
Nr pomiaru 1 2 3 4 5 6 7 8 9 10
IloĞü km (tys) 0 2 4 6 8 10 12 14 16 18
StopieĔ zuĪy-cia (%)
0 5 7 11 16 19 25 29 34 38
Po zastosowaniu algorytmu szarego modelu I rzĊdu otrzymano prognozĊ zuĪycia amortyzatora na najbliĪszych 10 tys. km przedstawioną w postaci wykresu
0 5 10 15 0 10 20 30 40 50 60 70 80 nr pomiaru w a rt o sc z u zycia Pomiary Prognoza
Rysunek 2. Prognoza zuĪycia amortyzatora na podstawie 10 obserwacji stanu Natomiast rzeczywiste wartoĞci zuĪycia przedstawiają siĊ nastĊpująco:
Tabela 3. Zaobserwowane wartoĞci zuĪycia amortyzatora w funkcji przebiegu pojazdu Nr pomiaru 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 IloĞü km (tys.) 0 2 4 6 8 10 12 14 18 20 22 24 26 28 30 WartoĞü zuĪy-cia amortyzato-ra 0 5 7 11 16 19 25 29 34 38 44 50 57 66 74
Po zestawieniu prognozy oraz wartoĞci otrzymanej z obserwacji otrzymujemy:
0.00 20.00 40.00 60.00 80.00 100.00 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 obserwacja - prognoza obserwacja
Rysunek 3. Zestawienie wartoĞci prognozowanych przy pomocy metody szarego modelu oraz wartoĞci zaobserwowanych
Z otrzymanego wykresu moĪna wnosiü, Īe:
– wartoĞci prognozowanego zuĪycia mają nieco wiĊkszą wartoĞü niĪ wartoĞci pochodzące z póĨ-niejszych obserwacji,
– róĪnice pomiĊdzy prognozą a wartoĞcią rzeczywistą nie przekraczają 7–8%,
– wartoĞü prognozowana odzwierciedla rzeczywisty stan zuĪycia amortyzatora, generując przy tym stosunkowo niewielki błąd,
– stosowanie metody szarego modelu (systemu) pierwszego rzĊdu moĪna uznaü za właĞciwe i uza-sadnione [4]
3.3. Zastosowanie metody szarego systemu pierwszego rzĊdu do prognozowania głĊbokoĞci bieĪnika w funkcji przebiegu samochodu osobowego [3]
Kolejne zastosowanie metoda szarego modelu pierwszego rzĊdu znajduje równieĪ na przykładzie nauk technicznych. Jako model przedstawiono prognozĊ zuĪycia bieĪnika opon, na podstawie przeprowadzonych obserwacji przy zastosowaniu algorytmu metody szarego systemu Deng’a. Pomiary zuĪycia bieĪnika przeprowadzano od momentu załoĪenia do samochodu osobowego, co 2 tys, km, do 20.tys. km.
Wyniki z przeprowadzonych pomiarów przedstawiają siĊ nastĊpująco:
Tabela 4 Zaobserwowane wartoĞci zuĪycia bieĪnika opony samochodu osobowego w funkcji przejechanych kilometrów
Nr pomiaru 1 2 3 4 5 6 7 8 9 10 11 IloĞü przejechanych km (tys.) 0 2 4 6 8 10 12 14 16 18 20 WartoĞü zuĪycia bieĪnika [mm] 0 0,21 0,38 0,56 0,76 0,99 1,31 1,59 1,86 2,21 2,58
Stosując algorytm szarego modelu (systemu) pierwszego rzĊdu, celem uzyskania prognozy, otrzymujemy wykres, przedstawiający dane wejĞciowe wraz z wartoĞciami prognozowanymi wpostaci wykresu, wygenerowanego przez Ğrodowisko Matlab
0 2 4 6 8 10 12 14 16 0 1 2 3 4 5 6 7 8 nr pomiaru w a rt o sc zu zyci a b ie zn ika [ m m ] Pomiary Prognoza
Rysunek 4. Zestawienie prognozowanych wartoĞci zuĪycia bieĪnika w funkcji przejechanych kilometrów
Otrzymane wartoĞci z prognozy porównano z wartoĞciami zaobserwowanymi z bezpoĞrednich pomiarów. Wyniki zaprezentowano w poniĪszej tabeli͘
Tabela 5 Zestawienie wartoĞci prognozowanych oraz zaobserwowanych zuĪycia bieĪnika opony w funkcji przebiegu, badanie własne
Graficznym przedstawieniem danych z tabeli jest wykres. 0 2 4 6 8 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 numer pomiaru w a rt o sc zu zyci a b ie z n ik a opony [ m m ] obserwacja prognoza I
Rysunek 5. Zestawienie wartoĞci prognozowanych oraz zaobserwowanych zuĪycia bieĪnika opony w funkcji przebiegu
Jak widaü na wykresie 5, wartoĞci prognozowane nieco odbiegają od zaobserwowanych, jed-nak ich rozrzut jest stosunkowo niewielki, wynosi ok. 8%. Warto zauwaĪyü, Īe wraz z kolejną pro-gnozą generowany błąd prognozy ulega zwiĊkszeniu, zatem metoda ta powinna byü stosowana do prognozowania krótkookresowego, jeĪeli priorytetem jest jak najmniejszy błąd [3].
4. Prognozowanie stanu obiektu w diagnostyce wielouszkodzeniowej przy pomocy szarego systemu pierwszego rzĊdu
Na chwilĊ obecną wielowymiarowa diagnostyka elementów oraz całych maszyn jest na początku swojego rozwoju, stąd nie moĪna mówiü jeszcze o pewnych i wypracowanych, na drodze obserwacji, metodach [1]. Jednak w przypadku okreĞlania prognozy stanu maszyn, na podstawie przetwarzanych symptomów, zastosowanie znajdują sieci neuronowe, jako jedno z istotnych zagadnieĔ, związanych ze sztuczną inteligencją.
Szare systemy równieĪ mogą byü stosowane do diagnostyki wielowymiarowej. Z badaĔ wynika, Īe optymalna szerokoĞü ruchomego okna dla GM(1,1) jest rzĊdu w = 4–8, mimo, Īe np. błąd prognozy nie jest jeszcze minimalny [1]. Bowiem powyĪej tej szerokoĞci przedziału wpogarsza siĊ jakoĞü prognozy i wielkoĞci prognozowane mogą byü równe lub nawet mniejsze od obserwowanych, co w diagnostyce jest zupełnie nieprzydatne [1].ĝredni błąd prognozy w funkcji szerokoĞci Ğlizgającego siĊ okienka został przedstawiony na poniĪszym wykresie:
Wykres 6 WartoĞü błĊdu prognozy w funkcji szerokoĞci Ğlizgającego siĊ okienka dla wielouszko-dzeniowej diagnostyki maszyn
ħródło: Diagnostyka’2/2007(1). 5. Podsumowanie
Niniejsza publikacja wskazuje na szerokie zastosowanie metody szarego systemu (modelu) pierwszego rzĊdu. Na chwilĊ obecną, przy pomocy tej metody, dokonuje siĊ prognozowania na potrzeby diagnostyki zarówno jednowymiarowej, jak i wielowymiarowej
NaleĪy pamiĊtaü, Īe metoda Denga wymaga właĞciwej implementacji w Ğrodowisku Matlab izastosowania parametrów, które pozwolą na prognozĊ, która nie bĊdzie obarczona duĪym błĊdem, np. właĞciwej szerokoĞci okienka Ğlizgającego dla metody RW (rolling window) – zarówno zbyt duĪa, jak równieĪ zbyt mała jego szerokoĞü moĪe wpływaü niekorzystnie na dokładnoĞü prognozy [3]. Przyjmuje siĊ, Īe szerokoĞü okienka Ğlizgającego dla prognoz
wielowymiarowych w granicach 4–8 daje minimalny błąd prognozy [1], natomiast szerokoĞü okienka poza tymi granicami wpływa na tyle niekorzystnie na jakoĞü prognozy, Īe moĪe ona róĪniü siĊ z znacznym stopniu od stanu obserwowanego. Uzyskanie właĞciwej prognozy stanu pozwala na zaplanowanie dalszych działaĔ, nie tylko w zakresie eksploatacji (badanie stanu obiektów technicznych), lecz takĪe w przypadku działaĔ związanych z zagadnieniami ekonomicznymi, np. przy planowaniu wielkoĞci produkcji.
Jak wykazano, metoda szarego systemu pierwszego rzĊdu jest metodą uniwersalną, przy jej stosowaniu naleĪy pamiĊtaü o tym, by dane wejĞciowe miały charakter rosnący oraz monotoniczny. Istotnym zagadnieniem jest okreĞlenie optymalnej szerokoĞci okienka Ğlizgającego, celem uzyskania jak najmniejszego błĊdu prognozy, a zatem najwiĊkszej przydatnoĞci metody szarego systemu w diagnostyce maszyn. Jest to przedmiotem działaĔ w zakresie dalszego rozpoznawania metody szarego systemu.
Bibliografia
1. Cempel C., Tabaszewski M., „Zastosowanie teorii szarych systemów do modelowania i prognozowania w diagnostyce maszyn”, Diagnostyka’2/2007: s. 11–18.
2. Jurek K. „Elementy sztucznej inteligencji w badaniu stanu maszyn – teoria szarego systemu”, Bydgoszcz 2010.
3. Jurek K., Opis i zastosowanie metody szarego systemu do prognozowania stanu obiektu. Wykorzystanie tej metody zarówno na potrzeby nauk technicznych, jak równieĪ innych gałĊzi nauki. W pływ szerokoĞci Ğlizgającego siĊ okienka na dokładnoĞü prognozy”, Bydgoszcz 2010. 4. Jurek K. „Zastosowanie metody szarego systemu (modelu) pierwszego rzĊdu na przykładzie nauk technicznych oraz ekonomicznych. OkreĞlenie błedu generowanego przez tĊ metodĊ”, Bydgoszcz 2010.
5. Jurek K., ĩółtowski B., „A review of artifical intelligence machine condition with a wide-spread consideration of grey system theory”, Studies&Proceedings Polish Association for Knowledge Management, Bydgoszcz 2010, s. 80–87.
6. Mierzwiak R., Werner K., Pochmara A., „Zastosowanie teorii systemów szarych do prognozowania ekonomicznych szeregów czasowych”, IV Krakowska Konferencja Młodych uczonych, Sesja Plenarna Nauki Ekonomiczne, Kraków 2009, s. 627–632.
THE IMPLEMETATION OF THE GREY SYSTEMS THEORY UNDER GENERAL DIAGNOSIS
Summary
The consideration on the grey systems theory as a forecasting method – which evaluates the condition of particular object – can be applied not only among tech-nical science but also in the other non-techtech-nical disciplines such as social and natu-ral science.
This publication aims to assess the origin of the grey systems theory as well as its mathematical interpretation. At the basis of the previous studies, it is easily to state that the first order grey model is frequently used by the scholars and engineers due to its transparency and interdisciplinary nature.
Keywords: grey systems, vibration diagnostics, forecasting, modeling of condition
PracĊ zrealizowano w ramach projektu „Techniki wirtualne w badaniach sta-nu, zagroĪeĔ bezpieczeĔstwa i Ğrodowiska eksploatowanych maszyn”.
Numer projektu: WND-POIG.01.03.01-00-212/09. Krzysztof Jurek
Wydział InĪynierii Mechanicznej
Uniwersytet Technologiczno-Przyrodniczy w Bydgoszczy e-mail: krzys89@onet.eu