Elektrotehniški vestnik 78(3): 142-146, 2011 Existing separate English edition Zaznavanje napak in spremljanje čiščenja odpadnih voda na podlagi mehkega modela Dejan DovZan1, Vito Logar2, Nadja Hvala3, Igor Škrjanc4 1'2'4Faculty of Electrical Engineering, Tržaška 25, 1000 Ljubljana 3Inštitut Jožef Stefan, Jamova cesta 39, 1000 Ljubljana E-pošta: 1 dejan.dovzan@fe.uni-lj.si Povzetek. V clanku je predstavljen sistem za spremljanje in zaznavanje napak na procesu čistilne naprave. Zaradi nelinearnosti procesa smo kot osnovo za sistem spremljanja uporabili mehki model, s katerim lahko zelo dobro aproksimiramo nelinearnosti procesa in se s tem izognemo nepotrebnim laznim alarmom. Model je dobljen s pomočjo Gustafson-Kesselovega algoritma za rojenje. Odziv mehkega modela v normalnem obratovanju je primerjan s trenutnim dogajanjem v procesu. Na podlagi odstopanja odzivov lahko zaznamo napako. Proces v našem primeru simuliramo. Simulirane pa so tudi napake na senzorjih. Merjeni signali so: vhodna koncentracija amonijaka in koncentracija kisika v prvem aerobnem tanku ter temperatura, koncentracija kisika in amonijaka v drugem tanku. Prikazani so rezultati spremljanja procesa in odkrivanja napak na podlagi mehkega modela. Ključne besede: mehko rojenje, graditev mehkega modela, cišcenje odpadnih voda, spremljanje procesa, zaznavanje napak Monitoring and sensor fault detection in a waste-water treatment process based on a fuzzy model In this paper, monitoring and sensor fault detection in a wastewater treatment process are discussed. Monitoring is based on the Takagi-Sugeno fuzzy model of a plant process obtained by using Gustafson-Kessel fuzzy clustering algorithm.The paper also explains the principle of the Takagi-Sugeno fuzzy model. The main idea is to cope with the non-linearity of a monitored process. The output of the fuzzy-model in normal operation regime is compared with the current behavior. If the fault-detection index exceeds a certain predefined value (the fault-tolerant index), an alarm is triggered. The data treated in this paper are obtained with a simulation model of a waste-water treatment plant and by simulating sensor faults. The signals to be measured in the process monitoring are the following: influent ammonia concentration, dissolved-oxygen concentration in the first aerobic reactor tank, temperature, dissolved-oxygen concentration and ammonia concentration in the second aerobic reactor. 1 UVOD Spremljanje procesov skupaj z zaznavanjem napak in diagnozo je v zadnjem casu postalo eno popularnejsih področij procesne avtomatike. Za spremljanje procesa in odkrivanje napak se pogosto uporabljajo tehnike na podlagi modela procesa, ekspertni sistemi in sistemi za razpoznavanje vzorcev [1], [2]. Hitri senzorji in oprema za zajem podatkov nam omogocajo zajem velike kolicine podatkov, ki nam kazejo, kaj se dogaja s procesom. Za obdelavo in analizo zajetih podatkov je bilo razvitih veliko statisticnih metod [3], [4]. S pomočjo Prejet 22. junij, 2011 Odobren 25. julij, 2011 meritev lahko ves cas spremljamo procese in zaznavamo napake ter temu primerno prozimo alarme. V našem primeru obravnavamo zaznavanje napake na senzorju simuliranega procesa cišcenja odpadne vode [5], [6], [7]. Ti procesi so zaradi dnevnih, mesecnih in sezonskih vplivov izpostavljeni nenehnemu spreminjanju dinamike. Nanje vplivjo sprememba temperature, kolicina padavin in razlicšne obremenitve naprave. Teoreticšno modeliranje je pri cišcenju odpadnih voda zelo zapleten postopek, katerega rezultati so dvomljivi. Zato se za spremljanje takih procesov uporabljajo statisticšne metode, ki temeljijo na rudarjenju s podatki. (Cišcenje odpadnih voda je zelo nelinearen proces, zato je treba njegovo nelinearnost upoštevati v modelu. Za predprocesiranje podatkov smo v nasšem primeru uporabili algoritem za mehko rojenje podatkov. S tem smo se izognili napacnim alarmom, ki lahko nastanejo zaradi nelinearnosti procesa (ce bi uporabili linearni model za zaznavanje), saj z mehkim modelom lahko poljubno natancno aproksimiramo nelinearnosti. 1.1 Mehki model z Gustafson-Kesselsovim rojenjem V tem poglavju bomo opisali algoritme in metode, ki so bile uporabljene za analizo podatkov. Razlozen bo Gustafson-Kesselov algoritem rojenja in izpeljan ter identificiran Takagi-Sugeno mehki model procesa cišcenja odpadnih voda. 1.1.1 Gustafson-Kesselov algoritem mehkega rojenja: Gustafson-Kesselov algoritem rojenja nam omogocša, da definiramo roje razlicšnih oblik, kar nam pri takih procesih, kot je obravnavani, pride prav. Matrika vhodnih podatkov je podana kot X G Rnxp. (1) Vektor vhodnih podatkov v trenutku k je definiran kot Xfc = [xfci,... ,xfcpj, xfc G Rp. Nabor n meritev označimo kot X = {xk | k = 1, 2,... ,n} in je predstavljen kot matrika n x p: X= xii X12 X21 X22 x1p X2p (2) (3) (4) Xn1 xn2 . . . xnp Čilj rojenja je razdeliti nabor podatkov X v c pod-mnoZic, ki jih imenujemo roji. Mehko razdeljen nabor podatkov X je skupek mehkih podmnozic (rojev) {A, | 1 < i < c}. Roji so definirani s pripadnostnimi funkcijami, ki so implicitno vsebovane v pripadnostni matriki U = [^jfc] G Rcxn. i-ta vrsta matrike U vsebuje vrednosti pripadnostne funkcije i-tega roja A, iz množice podatkov X. Matrika U izpolnjuje naslednje pogoje: pripadnostne stopnje so realna števila iz intervala od nic do ena G [0,1], 1 < i < c, 1 < k < n;), pripadnost vzorca xk k uniji rojev je ena (^°=1 = 1, 1 < k < n;), noben roj ni prazen oziroma ne vsebuje vseh podatkov (0 < J2n=1 Mifc < n, 1 < i < c). Zgoraj nasšteto pomeni, da matrika U pripada mnozšici, ki je definirana kot: M = {U G Rcxn | Mifc G [0,1], Vi, k; c n (5) = 1,Vk; 0 < ^^ifc < n, Vi}. i=1 k = 1 V nasšem primeru je matrika pripadnostni (U) dobljena z mehkim "c-means" algoritmom, ki temelji na Ma-halanobisovi normi. Algoritem dobimo z minimizacijo kriterijske funkcije ob upoštevanju omejitev iz En. 5: c n c n ,m J 2 J(X, U, V, A) = £ ]T .m^ + A £ ]T (Mik -1), i=1 k=1 i=1 k=1 (6) kjer je U matrika pripadnostni mnozice podatkov X, V je vektor centrov rojev V = [v1, v2,..., vc], v G Rp, (7) Ci En=1 Mik (Xfc- v»)(Xfc- v»)T En .m k = 1 Mifc Slednje nam omogoca detekcijo hiperelipsoidnega roja v porazdelitvi podatkov. Čš e so podatki porazdeljeni ob nelinearni hiperpovrsšini, algoritem najde roje, ki so lokalni linearni aproksimatorji hiperpovrsšine. Prekrivanje rojev je definirano s parametrom m G [1, to). Za dolocšanje sštevila rojev lahko uporabimo funkcije za merjenje veljavnosti rojev ali pa jih iterativno dodajamo oziroma odstranjujemo glede na odstopanje modela. Parameter prekrivanja m vpliva na mehkost roja; od trdega (m = 1) do cisto mehkega (m ^ to). V nasšem primeru smo parameter nastavili m = 2, kar je tudi obicšajna vrednost tega parametra. 1.1.2 Koraki Gustafson-Kesselovega algoritma: Postopek Gustafson-Kesselovega algoritma za rojenje nabora podatkov X lahko opišemo z naslednjimi koraki: • Inicializacija Izberemo število rojev c, faktor prekrivanja m (v našem primeru m = 2) in napako, pri kateri se algoritem konca eend > 0 (v našem primeru eend = 0.001). Inicializiramo matriko pripadnosti: U g M (nakljucno) in epoh r = 0. • Zanka r = r + 1 izračun centrov rojev: (r) V> = V-n ( (r)\m Efc=1 (vM(r ) Xk 1 < i < c. (8) Tn M(r) izračun kovariančne matrike in matrike A,: Cj = EL1 mZŽ (x& - vi)(Xfc - vi) En ,,m k=1 Mifc (9) A, = (pjdet (Cj))1/p C-1, 1 < i < c (10) izračun razdalje vzorcev do centrov rojev d2 "jfc (xfc - *,(r)) Aj (xfc - v(r)) , (11) 1 < i < c, 1 < k < n osvezitev matrike pripadnosti: (r) ce je djfc > 0, = f dik = M djk 2 (12) dk pa je kvadrat razdalje dk = (xfc - Vj)T A, (xfc - v,). Matrika Aj je definirana kot: A, = (pjdet (Cj))1/p C-1, kjer je p, = 1, i = 1,..., c in p enak številu merjenih spremenljivk. Cj je mehka kovariancšna matrika i-tega roja, definirana kot • dokler ni ||U(r) - U(r-1)|| < eend 1.1.3 Mehki model Takagi-Sugeno : Mehki TS-modeli aproksimirajo nelinearen sistem z interpolacijo lokalnih linearnih modelov. Vsak lokalni linearni model prispeva svoj delez k izhodu globalnega modela. Če imamo nabor vhodnih vektorjev iT X = [X1 ,X2, ...,Xn] in nabor pripadajocih izhodov Y = [y1,y2,... ,yn]T (13) m 1 potem je mehki model podan v obliki pravil R^ Ri : ce je xk iz Ai potem yk = ^i(xk), i = 1,..., c (15) Vektor xk je vektor vhodnih podatkov, je izhod lokalnega linearnega modela v trenutku k. Vektor xk pripada mehkim mnoZicam (Ai,...,Ac) z določeno pripadnostjo (xk) oziroma pik : R ^ [0,1]. Funkcije ) so lahko poljubne zvezne funkcije, čeprav so največkrat uporabljene linearne afine funkcije. Izhod globalnega modela izracunamo po enacbi: Vk Ec=i Mik&i(xk) (16) ^¿=1 pik En. (16) poenostavimo tako, da definiramo funkcijo: A(xk) = Pik i = 1, . . . , c (17) Vk =53 Pi (Xk )ti (Xk ), k =1,..., (18) Izhod lokalnega modela je obicajno definiran kot linearna kombinacija elementov podatkovnega vektorja: $i(xk) = Xk0i, i = 1,... ,c, = [°i1, °i(p+q)\ . (19) Vektor zmehcanih vhodov v trenutku k zapišemo kot: ^k = [A(xk )xk,.. .,Pc(xk)xk ], k = 1,...,n, (20) matriko zmehcšanih podatkov pa kot: = OT, i>T,..., ]. (21) Ce zapišemo matriko koeficientov celotnega nabora pravil kot eT = K ], (22) potem lahko izhod globalnega modela (En. (18)) zapišemo v matricni obliki: Vik = ^k ©. (23) Relacija izhoda modela in vhodnih podatkov za celoten nabor zapisšemo kot: Y = ^0, (24) kjer je Y vektor izhodov modela yk (k = 1,..., n) Y = [V 1,V2,...,Vn ] . Mehki model, podan z enacbo En. (23), imenujemo afini model Takagi-Sugeno. S tem modelom lahko s poljubno natancnostjo aproksimiramo poljubno funkcijo [8], [9], [10]. Splošnost modela lahko dokazemo s Stone-Weierstrassovim teoremom [11]. Ta pravi, da lahko vsako zvezno funkcijo aproksimiramo z razsširitvijo osnovnih mehkih funkcij [12]. 1.1.4 Ocena parametrov mehkega modela: Za oceno parametrov mehkega modela bo uporabljena metoda najmanjših kvadratov. Meritve zadošcajo nelinearni enacbi sistema: Vi = g(xi), i = 1,... ,n (26) Po Stone-Weierstrassevm teoremu obstaja za realno zvezno funkcijo g, definirano nad naborom Uc C Rp, mehki sistem f tako, da je Sc=1 pik Funkcija ^i(xk) nam podaja normirano pripadnost vzorca k posameznemu roju. S tem nam pove tudi, kaksšna je veljavnost pravila. Vsota funkcije po rojih je ena (EC=1 A(xk) = 1) ne glede na xk, ce je imenovalec ^i(xk) razlicen od nic. To pa je lahko zagotoviti s pravilno definicijo pripadnostnih funkcij nad roji. Z zdruzitvijo enacb (16) in (17) pridemo do slednje poenostavljene enacšbe za izhod globalnega modela: max |f (xi) - g(xi Xi£X < s, yi, (27) kjer je S > 0 poljubno majhna konstanta. Pri aproksima-ciji zveznih funkcij z mehkimi funkcijami iz razreda Fp, ki so definirane z En. (23), je treba opozoriti, da manjše vrednosti S povecajo število rojev c, da bi zadovoljili En. (27). Pri aproksimaciji funkcije z mehkim modelom je napaka med meritvami in vrednostmi izhoda mehkega modela definirana kot: Vi - f (Xi) = Vi - Vi, i = 1,...,n, (28) kjer Vi pomeni merjene izhode sistema, V k pa izhode mehkega modela v trenutku k. Parametre mehke funkcije (0) dobimo z minimizacijo vsote kvadratov napak: E- E' i=i (29) (25) = (Y - Y)T(Y - Y) == (Y - ^0)T(Y - ^0). Parameter 0 dobimo kot parcialni odvod §§ = 0: 0 = (^T Y. 2 Biološka (Čistilna naprava za ČIŠČENJE ODPADNIH VODA Procesi cišcenja odpadnih voda so veliki nelinearni sistemi, ki so izpostavljeni velikim motnjam predvsem v toku in obremenitvi. Hkrati pa se spreminja tudi sestava odpadne vode, ki pride v cistilno napravo. Za nepristransko primerjavo razlicšnih sistemov oziroma strategij vodenja in sistemov spremljanja je bil razvit enoten simu-lacijski model. Sestoji iz petih zaporedno vezanih reaktorjev, ki se nato iztekajo v desetplastni sekundarni tank, kjer se tvorijo usedline. Model procesa, pripadajoce enacbe in podroben opis najdemo na spletni strani http : //www.ensic.inpl - nancv.fr/COSTWWTP/. V nasšem primeru smo spremljali del sistema, kjer je odpa-dana voda cišcena (tanki). Shema procesa je prikazana na sliki 1. e n 2 tanki anoxic tanks T T T T T T T T T T T T aeration tanks 0 0.5 1 1.5 2 Sample H O 5 O 5 0 2.5 3 0 0.5 1 1.5 2 Sample 2.5 3 4 Slika 2: Čeloten nabor meritev. Koncentracija amonijaka v pritoku v sistem CNH4Nin, koncentracija kisika v prvem tanku C O 2 ter koncentracija kisika CO 2 in amonijaka CNH4N out v drugem tanku. Slika 1: Shema simuliranega procesa. Za izračun mehkega modela obravnavanega dela procesa smo uporabili naslednje spremenljivke: koncentracija amonijaka v vhodnem toku v proces Qin, ki je oznacena kot CNH4N, koncentracija stopljenega kisika v prvem tanku Cq , koncentracija stopljenega kisika v drugem tanku Cq2 in koncentracija amonijaka v drugem tanku CNH4Nout. Mehki model je bil zgrajen tako, da aproksimira odvisost koncentracije amonijaka v drugem tanku od preostalih merjenih spremenljivk: CNH4NoUt (k) = G (CnH4N„ (k),CQ2 (k),C%2 (k)) , (30) kjer G oznacuje nelinearno povezavo med merjenimi spremenljivkami. Prvih 15000 vzorcev (s casom vzorcenja Ts = 120s) je bilo uporabljenih za gradnjo (identifikacijo) mehkega modela Takagi-Sugeno. Ob 17000. vzorcu zacne pocasi narašcati napaka na senzorju, ki je odpravljena ob 18000. vzorcu. Senzor za merjenje koncentracije amonijaka v drugem tanku CNH4N out seje torej pokvaril. K nominalnemu signalu senzorja smo dodali eksponentni signal za simulacijo napake. Čeloten set meritev je prikazan na sliki 2. 100 r 50 - Sample Slika 3: Verifikacija mehkega modela, kjer je CNH4Nout izhod modela, CNH4Nout pa izhod procesa. napake je definiran kot: C NH4N o -c. NH 4 N o C (31) NH4N o Ko indeks detekcije preseže indeks tolerance napake, se sproži alarm. Indeks tolerance je defeniran kot največji indeks napake, pomnožen s konstanto: ftoi = Y max f. Največji indeks napake smo zaznali pri učenju modela na podatkih, kjer ni bilo simulirane napake. V našem primeru smo pomnožili maksimalno vrednost indeksa napake s konstanto y = 1,5. To pomeni, da smo za tolerančni indeks dobili vrednost ftol = 0,15. Slika 4 prikazuje tolerančni indeks, indeks detekcije napake in signal za alarm. Napako, ki nastane pri vzorcu 17000, je predlagani sistem za spremljanje procesa detektiral pri 17556. vzorcu. Detekcija napake je zakasnjena, kar je pri napakah, ki naraščajo počasi običajno. Ji Sample 1.5 Sample Mehki model smo zgradili na prvih 15000 vzorcih. Izhod mehkega modela CNH4Nout in izhod procesa CNH4Nout sta prikazana na sliki 3. Indeks za detekcijo Slika 4: Indeks detekcije napake, tolerancni indeks ftoi dejanska in detektirana napaka. q 4 2 0 5000 10000 15000 2 0 10 10 2 x 10 x 10 3 Sklep V clanku smo obravnavali detekcijo napake na senzorju in spremljanje čiščenja odpadnih voda. Sistem za detekcijo je bil realiziran z mehkim modelom Takagi-Sugeno. Model smo identificirali s pomočjo Gustafson-Kesslovega algorita rojenja, parametri pa so bili dobljeni s pomočjo najmanjših kvadratov. Razvit sistem je bil testiran na modelu procesa, kjer je bila simulirana napaka na enem od senzorjev. Za gradnjo mehkega modela smo uporabili naslednje meritve: vhodno koncentracijo amo-nijaka in koncentracijo kisika v prvem tanku ter temperaturo, koncentracijo kisika in koncentracijo amonijaka v drugem. Zaradi narave napake, simulirane na senzorju koncentracije amonijaka na drugem tanku, je bila le-ta detektirana brez lazšnih alarmov in z majhno zakasnitvijo. Glede na to, da je dinamika sistema izpostavljena spreminjanju zaradi razlicšnih obremenitev, padavin itd., bomo poskusšali sistem za detekcijo napake implementirati z rekruzivno mehko identifikacijo, ki bo sposobna prilagajati dinamiko sistema tem spremambam, hkrati pa bo sistem sposoben detektirati tudi napake. Vito Logar je diplomiral leta 2004 in doktoriral leta 2009 na Fakulteti za elektrotehniko v Ljubljani. Podrocje njegovega raziskovanja zajema napredno analizo EKG signalov in modeliranje industrijskih procesov. Trenutno dela na procesu elektricne oblocne peci in optimizaciji stroškov njenega delovanja. Nadja Hvala je diplomirala leta 1985, magistrirala leta 1988 in doktorirala leta 1992, vse na Fakulteti za elektrotehniko Univerze v Ljubljani. Kot znanstvena sodelavka je zaposlena na Institutu Jozef Stefan v Ljubljani. Ukvarja se z modeliranjem, simulacijo, optimizacijo in nacrtovanjem vodenja kemijskih in biokemijskih procesov. Igor Skrjanc je diplomiral leta 1988, magistriral leta 1991 in doktoriral leta 1996 na Fakulteti za elektrotehnik v Ljubljani. Njegovo raziskovalno podrocje obsega adaptivne, mehke in mehke adaptivne regulacijske sisteme. Leta 2007 je prejel Vodovnikovo nagrado, leta 2008 Zoisovo nagrado Republike Slovenije za dosežke na raziskovalnem področju pametnega vodenja. Za obdobje med 2009 in 2011 je prejel Humboldtovo štipendijo za izkušene raziskovalce za raziskovalno delo na Univerzi v Siegenu. Literatura [1] Chen, J., Liao, C. M., "Dynamic process fault monitoring based on neural network and PCA". Jour. of Process Control, vol. 12, pp. 277-289, 2002. [2] Klancar, G., D. Juricic, R. Karba., "Robust fault detection based on compensation of the modelling error". International journal of Systems Science, vol. 33, no. 2, pp. 97-105, 2002. [3] Johnson, R. A., Wichern, D. W., Applied Multivariate Statistical Analysis, Prentice-Hall, New Jersey, 1992. [4] Daszykowski, M., Walczak, B., Massart, D. L., "Projection methods in chemistry". Chemometrics and Intelligent Laboratory Systems, vol. 65, pp. 97-112, 2003. [5] Vrecko, D., Hvala, N., Kocijan, J., Zec, M. "System analysis for optimal control of a wastewater treatment benchmark". Water sci. technol., vol. 43, pp. 199-206, 2001. [6] Vrecko, D., Hvala, N., Kocijan, J., "Wastewater treatment benchmark : What can be achieved with simple control?". Water sci. technol., vol. 45, pp. 127-134, 2002. [7] Hvala, N., Vrecko, D., Burica, O., StraZar, M., Levstek, M., "Simulation study supporting wastewater treatment plant upgrading". Water sci. technol., vol. 46, pp. 325-332, 2002. [8] Kosko, B., "Fuzzy Systems as Universal Approximators", IEEE Transactions on Computers, vol. 43, no. 11, pp. 1329-1333, 1994. [9] Ying, H. GC., "Necessary conditions for some typical fuzzy systems as universal approximators", Automatica, vol. 33, pp. 1333-1338, 1997. [10] Wang, L.-X., Mendel, J. M., "Fuzzy basis functions, universal approximation, and orthogonal least-squares learning", IEEE Trans. Neural Networks, vol. 3, no. 5, pp. 807-814, 1992. [11] Goldberg, R. R., Methods of Real Analysis, John Wiley and Sons, 1976. [12] Lin, C-H., "Siso nonlinear system identification using a fuzzy-neural hybrid system", Int. Jour. of Neural Systems, vol. 8, no. 3, pp. 325-337, 1997. Dejan Dovzan je diplomiral leta 2008 na Fakulteti za elektrotehniko. Trenutno je zaposlen kot mladi raziskovalec v laboratoriju Laboratorij za modeliranje simulacijo in vodenje na Fakulteti za elektrotehniko. Njegovo raziskovalno podrocje obsega adaptivno mehko vodenje, pre-diktivno vodenje, mehke modele, detekcijo napak s pomocšjo mehkih modelov in rekurzivno mehko identifikacijo.