ISSN 2350-6113 P o mu r s k a L e t n i k 1, l e t o 2014, št. 2 Obzorja Izdajatelj Pomurska akademsko znanstvena unija Diana Gregor ŽIVLJENJSKI KROG RAČUNALNIK, Domen LiDAR Svetec, KI Mongus: - Tanja SIMULIRANJE Silva PAPIRJA JE König, IN PAPIRNE SPREMENIL ALGORITMI Bagar: Janja SVET NAD - EMBALAŽE Borut VELIKIMI ALTERNATIVNA ELEKTROENERGETSKEGA Zule, Gorazd - Jože Žalik, Niko KOLIČINAMI GORIVA SISTEMA – - Nemec: Lukač, PODATKOV Rafael ZAKAJ Golob: IN Mihalič: KAKO? - Zoran Žunič: MEŠANA METODA ROBNIH IN KONČNIH ELEMENTOV ZA REŠEVANJE HITROSTNO VRTINČNE FORMULACIJE NAVIER-STOKESOVIH ENAČB - Iztok Fister, Iztok Fister, Janez Brest, Iztok Fister Jr., Karin Ljubič: Jeke, ALGORITMI Matjaž KAPLJEVITEGA PO VZORIH Hriberšek: FILMA - IZ NARAVE NUMERIČNA Matej Zadravec, - Matej Zadravec, SIMULACIJA Matjaž Uroš UPARJANJA Hriberšek, Jure Ravnik: LOČEVANJE MAGNETNIH DELCEV IZ TEKOČINE - Milan Šernek, Maksimilijan Mravljak: ADHEZIJA PRI LEPLJENJU LESA Z JEKLOM P O M U R S K A A K A D E M S K O Z N A N S T V E N A U N I J A Pomurska obzorja Letnik 1, leto 2014, številka 2 ISSN 2350-6113 Murska Sobota, 2014 Kazalo IMPRESUM UVODNIK 47 47 TEHNIKA ■ ŽIVLJENJSKI KROG PAPIRJA IN PAPIRNE EMBALAŽE Diana Gregor Svetec et al. 48 ■ RAČUNALNIK, KI JE SPREMENIL SVET Jože Nemec 51 ■ ALGORITMI NAD VELIKIMI KOLIČINAMI PODATKOV LiDAR Borut Žalik, Niko Lukač, Domen Mongus 54 ■ ALTERNATIVNA GORIVA Tanja Bagar 63 ■ SIMULIRANJE ELEKTROENERGETSKEGA SISTEMA – ZAKAJ IN KAKO Rafael Mihalič 66 ■ MEŠANA METODA ROBNIH IN KONČNIH ELEMENTOV ZA REŠEVANJE HITROSTNO VRTINČNE FORMULACIJE NAVIER-STOKESOVIH ENAČB Zoran Žunič 72 ■ ALGORITMI PO VZORIH IZ NARAVE Iztok Fister et al. 76 ■ NUMERIČNA SIMULACIJA UPARJANJA KAPLJEVITEGA FILMA Matej Zadravec, Uroš Jeke, Matjaž Hriberšek 82 ■ LOČEVANJE MAGNETNIH DELCEV IZ TEKOČINE Matej Zadravec, Matjaž Hriberšek, Jure Ravnik 87 BIOTEHNIKA ■ ADHEZIJA PRI LEPLJENJU LESA Z JEKLOM Milan Šernek, Maksimilijan Mravljak 92 Uvodnik Impresum ZDRAVJE V ŠIRŠEM POMENU Internet: http://www.pazu.si e-mail: pazu@pazu.si ISSN 2350-6113 Naslov publikacije: POMURSKA OBZORJA Letnik 1 Leto 2014 Številka 2 Izdajatelj: Združenje Pomurska akademsko znanstvena unija Uredništvo: Odgovorni urednik pom. akad. dr. Mitja Slavinec Glavni urednik pom. akad. dr. Milan Svetec Tehnični urednik Zoran Wolf Uredniški svet: akad. pom. akad. dr. Anton Vratuša pom. akad. dr. Damir Josipovič pom. akad. dr. Albina Nećak Lük pom. akad. dr. Vesna Kondrič Horvat pom. akad. dr. Darja Senčur-Peček pom. akad. dr. Mitja Lainščak pom. akad. dr. Mirjam Sepesy Maučec pom. akad. dr. Rafael Mihalič akad. pom. akad. dr. Igor Emri Oblikovanje naslovnice: doc. mag. Tilen Žbona Tisk: Tiskarna aiP Praprotnik d.o.o. Tavčarjeva ulica 14, Černelavci Naslov izdajatelja in uredništva: PAZU – Uredništvo revije Pomurska obzorja Lendavska ulica 5a, 9000 Murska Sobota Datum natisa: maj 2014 Naklada: 200 izvodov V drugi številki Pomurskih obzorij objavljamo preostale članke z 11. znanstvene konference Pomurske akademije PAZU, ki se je med drugim ponašala z zelo pomenljivim dejstvom. Na konferenci je sodelovala preko četrtina vseh članic in članov PAZU. S tem je Konferenca z velikim naskokom potrdila, da ima vlogo najpomembnejšega znanstveno raziskovalnega dogodka v Pomurju. To je le še dodatni razlog, zakaj smo naše osrednje stičišče za razprave strokovnjakov in znanstvenikov, vsebinsko nadgradili tudi z revijo, katere drugo številko imate v roki. V prejšnji številki smo predstavili članke, neposredno povezane s plenumom 11. Konference PAZU Zdravo Pomurje za Zdrav razvoj. Naš motiv in ambicije z izbiro rdeče niti konference pa so bili veliko širši. Zdravja nismo hoteli obravnavati zgolj kot odsotnost bolezni pri ljudeh, temveč ga obravnavati na čim več različnih področjih, kot so zdravo okolje, zdrava predelava hrane, ki je v Pomurju pogosto kar samo po sebi umevna, zato to našo prednost pogosto kar spregledamo. S plenumom smo na osnovi širokega znanstvenega potenciala, ki ga PAZU zagotavlja, želeli opozoriti na čim več možnosti, ki nam jih nudi naša pokrajina. V drugi številki predstavljamo prispevke z 11. PAZU konference, ki se nanašajo na področje naravoslovja in tehnike. Tako je pomurska akademikinja dr. Diana Gregor Svetec s sodelavci, predstavila možnosti večje skrbi za okolje z recikliranjem papirja, dr. Jože Nemec nas je spomnil začetkov osebnih računalnikov, dr. Borut Žalik je s sodelavci, prav na primeru iz Pomurja, predstavil algoritme, ki operirajo z velikim številom podatkov. Tudi prispevek dr. Tanje Bagar o alternativnih gorivih je povezan z zdravjem našega okolja, dr. Rafael Mihalič pa nas pouči, kako funkcionira sistem, brez katerega si življenja več ne moremo predstavljati, to je elektroenergetika. Dr. Zoran Žunič predstavlja posebno metodo za reševanje problemov s področja hidrodinamike, dr. Iztok Fister in sodelavci, pa so se pri pripravi algoritmov zgledovali po naravi. Dr. Matej Zadravec s sodelavci piše o neslutenih možnostih numeričnih simulacij, revijo pa zaključuje članek rednega gosta na konferencah, dr. Milana Šerneka o možnostih lepljenje lesa z jeklom, kar lesu izredno poveča uporabnost. pom. akad. dr. Mitja Slavinec predsednik PAZU Revija Pomurska obzorja izhaja dvakrat letno. Revija je brezplačna. Pomurska obzorja 1/ 2014/ 2 | 47 Tehnika Diana Gregor Svetec1, Silva König1, Janja Zule2, Gorazd Golob1 Življenjski krog papirja in papirne embalaže POVZETEK Evropska uredba o embalaži in odpadni embalaži določa ukrepe, katerih poglavitni cilj je preprečiti nastajanje odpadne embalaže. Druga osnovna načela uredbe pa so ponovna uporaba embalaže, recikliranje in druge oblike predelave odpadne embalaže in s tem zmanjšanje dokončnega odstranjevanja odpadkov. Uredba o odpadkih zahteva, da se spodbuja visokokakovostno recikliranje in sistemi ločenega zbiranja odpadkov, primerni za doseganje potrebnih standardov kakovosti recikliranja, če je to tehnično, okoljsko in ekonomsko izvedljivo. Teoretična zgornja meja pri zbiranju in predelavi papirja ter papirne embalaže primerne za recikliranje je 81 %. Cilj Evropske skupnosti je doseči stopnjo recikliranja višjo od 70 % in uporabiti 50 % delež recikliranih vlaken v papirju in kartonu. Eden izmed ciljev projekta “Eko(loški) krogotok papirja” je povečati ozaveščenost in izdelati priporočila za optimalni postopek zbiranja papirja in papirnih izdelkov po posameznih regijah ter pripraviti izhodišča za novo transnacionalno strategijo vezano na uporabo recikliranih papirjev in papirne embalaže. Ključne besede: papir, papirna embalaža, življenjski krog, CE projekt. Projekt EcoPaperLoop Projekt EcoPaperLoop je sofinanciran s strani Evropske Unije/Evropskega regionalnega razvojnega sklada (ERDF) in lokalnih projektnih partnerjev. Vanj je vključenih 5 držav oz. 10 partnerjev iz osrednje Evrope, in sicer iz Slovenije dva, Italije in Nemčije po trije, iz Madžarske in Poljske pa po eden (1). Glavni namen projekta je ozaveščati javnost o pomenu celotnega ekološkega krogotoka papirja in papirne embalaže ter zagotoviti orodja za učinkovito zbiranje, ponovno uporabo oz. recikliranje. V regijah Srednje Evrope predstavlja recikliran papir pomemben surovinski vir. Stopnja recikliranja je kljub vsem trudom ozaveščanja javnosti še vedno izjemno nehomogena in se velikokrat izvaja na drugi lokaciji, kot je bil papir proizveden. Pri tem sta ključnega pomena ekološko oblikovanje in sistematično zbiranje, ki morata biti osnovana na podlagi mednarodnega vzajemnega sodelovanja vseh članic Srednje Evrope (2). Življenjski krog papirja Najpomembnejša surovina za proizvodnjo papirja so vlaknine. V svetu pridobivajo kar 90 % vlaknin za papir iz različnih vrst lesa iglavcev (smreka, jelka, bor) in listavcev (bukev, breza, topol). Glede na sestavo gozda in razpoložljivost posameznih vrst lesa uporabljajo v različnih državah tudi nekoliko različne deleže lesnih vrst za pridobivanje vlaknin. Kljub temu je povprečni delež lesa iglavcev precej večji od 1 48 | Univerza v Ljubljani, Naravoslovnotehniška fakulteta, Oddelek za tekstilstvo, Snežniška 5, Ljubljana 2 Inštitut za celulozo in papir, Bogišićeva 8, Ljubljana Pomurska obzorja 1/ 2014/ 2 deleža lesa listavcev. Les iglavcev in listavcev se razlikuje po gostoti (gostota lesa listavcev je večja kot iglavcev), dolžini vlaken (iglavci imajo v povprečju daljša vlakna) in debelini stene (iglavci imajo daljše stene) (3, 4). Sveže vlaknine prihajajo v papirnico običajno v suhem stanju v obliki velikih skladovnic ali bal večjega formata. Iz skladišča potujejo te bale po tekočem traku v razpuščevalnik za pripravo papirovine ali papirne snovi. Razpuščevalnik je velika posoda z mešalom na dnu, v katero natočijo vodo in ji dodajo določeno število bal celulozne ali druge vlaknine. Med intenzivnim mešanjem se vlaknine hitro razvlaknijo v vlakninsko kašo. Zaradi prisotnosti nerazpuščenih kosmov vlaken v razpuščeni snovi vodijo papirno snov še skozi razkosmovalnik, v katerem se kosmi zaradi velike hitrosti snovi med noži razkosmijo. Najvažnejši postopek pri pripravi snovi je mletje, saj pri tem procesu dobijo vlaknine želeno dolžino in strukturo. Včasih so mleli vlaknine v holandcu, danes pa v rafinerju, ki ima to prednost, da postopek poteka kontinuirno. Zmleta vlakna in dodatke (polnila, barvila) dodajajo s pomočjo dozirne naprave v mešalno kad. Papirna snov se mora ves čas mešati, da se njena konsistenca ne spremeni. Iz mešalne kadi prečrpavajo papirno snov v strojno kad od tam pa na papirni stroj (3, 4). Da bi prihranili surovine in energijo, lahko papir tudi recikliramo. To pomeni, da papirna vlakna ponovno uporabimo, namesto da bi primarna vlakna pridobili iz lesa (5). V papirnicah papir in karton pripravijo za obdelavo – odstranijo lepilo, vezavo, kovinske dele in smeti ter surovine primerno predelajo (6, 7). Star papir v obliki bal ali v razsutem stanju običajno transportirajo na tekočem traku v razpuščevalnik z vročo vodo ali vodo iz krogotoka papirnega stroja, kjer se dobro premeša. Večkrat dodajajo kemikalije, kot so dodatki za razsivitev (3, 8). Diana GREGOR SVETEC et al.: ŽIVLJENJSKI KROG PAPIRJA IN PAPIRNE EMBALAŽE Razsivitveni postopek je kemijsko-mehanski postopek za odstranjevanje tiskarskih barv s potiskane površine starega papirja. Postopek je v osnovi sestavljen iz treh glavnih faz: odstranjevanje tiskarske barve z vlaken s pomočjo kemikalij, odstranjevanje tiskarske barve s snovi s pomočjo flotacije ali izpiranja ter dokončna priprava recikliranih vlaken z zgoščevanjem, nevtralizacijo in z morebitnim beljenjem (3, 7). Sledi ponovna izdelava papirja, tokrat iz recikliranih vlaken, po konvencionalnem postopku. Poraba in proizvodnja papirja v Evropi Poraba papirja se je skozi leta močno povečala. Povprečna letna poraba v Evropi znaša 220 kg na prebivalca, pri čemer je v Nemčiji največja z 254 kg. V Sloveniji je leta 2010 znašala 185 kg. Največji porabniki papirja so Američani z 288 kg, manj pa so doslej porabili papirja na Kitajskem s 55 kg na prebivalca letno (9). Evropska konfederacija proizvajalcev papirja (CEPI) predstavlja 95 % proizvodnje evropske celulozne in papirne industrije. Med članicami je tudi Slovenija. Za izdelavo papirja, kartona oz. lepenke se je v letu 2012 v državah CEPI porabilo 41,9 milijona ton (40,0 %) lesnih vlaken, 0,4 milijona ton (0,4 %) drugih vlaken, 46,8 milijona ton (44,7 %) papirja za recikliranje in 15,6 milijona ton (14,9 %) nevlakenskega materiala (10). Skozi leta se je delež papirja za recikliranje, gledano na celoten porabljen papir v Evropski uniji ter Norveški in Švici, povečeval. Leta 2012 je delež znašal 71,7 %. Od celotne količine papirja za recikliranje je 66 % namenjenega embalažnemu papirju, 17 % časopisnemu, 7 % grafičnemu (mehanski in brezlesni), 7 % gospodinjskemu in sanitarnemu ter 3 % preostalemu papirju. Časopisni in embalažni papir sta večinoma sestavljena iz papirja za recikliranje (91 oz. 76 %), gospodinjski in sanitarni papir prav tako vsebujeta veliko količino papirja za recikliranje (46 %), medtem ko se za grafični papir še vedno večinoma uporablja sveža vlakna (samo 11 % papirja za recikliranje). Papir za recikliranje je večinoma sestavljen iz valovitega kartona in kraft papirja (45 %), sledijo časopisi in revije (25 %), papirji mešanih razredov (20 %) in papirji višjega razreda (10 %) (10). Zakonodaja V Evropi in Sloveniji so sprejete naslednje pomembnejše direktive in uredbe na področju krogotoka proizvodnje, uporabe, zbiranja in recikliranja papirja: Direktiva o odpadkih, Uredba o odpadkih, Direktiva o embalaži in odpadni embalaži in Uredba o zelenem javnem naročanju. Slovenija je z vstopom v Evropsko Unijo svojo okoljsko zakonodajo postopoma uskladila z evropsko, ki v Uredbi o odpadkih določa naslednjo hierarhijo pri ravnanju z odpadki: 1. preprečevanje, 2. priprava za ponovno uporabo, 3. recikliranje, 4. drugi postopki predelave (npr. energetska predelava) in 5. odstranjevanje odpadkov. Uredba določa, da je potrebno zagotoviti ukrepe za spodbujanje visokokakovostnega recikliranja in sistemov ločenega zbiranja odpadkov, primernih za doseganje potrebnih standardov kakovosti recikliranja, če je to tehnično, okoljsko in ekonomsko izvedljivo. Odpadke iz papirja, kovine, plastike in stekla je treba zbirati ločeno. Eden izmed ciljev uredbe je, da se do leta 2020 priprava za ponovno uporabo in recikliranje odpadnega papirja, kovin, plastike in stekla iz gospodinjstev, v čim večji meri pa tudi iz drugih virov, kjer gre za tokove odpadkov, podobne odpadkom iz gospodinjstev, povečata na najmanj 50 odstotkov skupne teže (11). Cilj Direktive o embalaži in odpadni embalaži (94/62/ES) je uskladitev nacionalnih predpisov za ravnanje z embalažo in odpadno embalažo, da se, po eni strani, prepreči njen kakršenkoli vpliv na okolje vseh držav članic kot tudi tretjih držav ali da se tak vpliv zmanjša, s čimer se zagotovi visoka raven varstva okolja, ter, po drugi strani, da se zagotovi delovanje notranjega trga in se preprečijo trgovinske ovire kot tudi izkrivljanje in omejevanje konkurence v Skupnosti. V ta namen direktiva določa ukrepe, katerih poglavitni cilj je preprečiti nastajanje odpadne embalaže, druga osnovna načela pa so ponovna uporaba embalaže, recikliranje in druge oblike predelave odpadne embalaže, in s tem zmanjšanje dokončnega odstranjevanja takšnih odpadkov (12). Namen Uredbe o zelenem javnem naročanju je zmanjšati negativen vpliv na okolje z javnim naročanjem okoljsko manj obremenjujočega blaga, storitev in gradenj in dajanje zgleda zasebnemu sektorju ter potrošnikom. Uredba določa okoljske zahteve tudi za pisarniški papir in higienske papirnate proizvode (13). Zakonodaja vsebuje tudi določena neskladja. Tako na primer v zgoraj omenjenih zakonih recikliranje materialov spodbujajo, določeni standardi (EN 643) že uporabljen papir priznavajo kot dragoceno surovino, medtem ko v določeni zakonodaji temu ni tako. Direktiva 2001/77/ES na primer spodbuja proizvodnjo električne energije iz obnovljivih virov energije (iz biomase, vključno z materialom iz gozdov in sorodnih industrij), s čimer opredeljuje celulozne in papirne ostanke kot obnovljiv vir energije (14). Direktiva o odpadkih (2008/98/ES) predstavlja možnost, da nekateri določeni odpadki prenehajo biti odpadki, ko so predelani, vključno z recikliranjem, in izpolnjujejo določena merila (15). JRC-IPTS (Joint Research Centre – Institute for Prospective Technological Studies) je objavil tehnični predlog o pogojih za prenehanje statusa odpadka za odpadni papir skladno z omenjeno direktivo o odpadkih. Tehnični predlog zajema analize recikliranja odpadnega papirja ter potencialnih ekonomskih, okoljskih in pravnih posledic konca statusa odpadka za odpadni papir (16). Evropska deklaracija o recikliranju papirja narekuje, da mora zbiranje papirja ostati vsaj pri trenutni visoki stopnji v državah, kjer že dosega stopnjo nad 70 % in mora narasti v državah, kjer je stopnja nižja od 60 %. Cilj evropske deklaracije o recikliranju papirja je doseči 70 % stopnjo recikliranja do leta 2015. Pri tem je potrebno vedeti, da je teoretična meja stopnje recikliranja papirja in kartona okoli 81 %. Medtem ko 19 % papirnih izdelkov zaradi tehničnih razlogov ni možno zbrati Pomurska obzorja 1/ 2014/ 2 | 49 Diana GREGOR SVETEC et al.: ŽIVLJENJSKI KROG PAPIRJA IN PAPIRNE EMBALAŽE (npr. arhivi knjižnic) ali reciklirati (npr. higienski papirji, kavni filtri) (17). Literatura 1. MOŽINA, K., GREGOR-SVETEC, D. in KÖNIG, S. Seminar recikliranja : Darmstadt. Grafičar, Ljubljana, 2013, št. 2, str.20. 2. KÖNIG, S., BLAZNIK, B., GOLOB, G., MOŽINA, K., URBAS, R., VRABIČ BRODNJAK, U. in GREGORSVETEC, D. 2. dogodek EcoPaperLoop : Ljubljana. Grafičar, Ljubljana, 2013, št. 4, str. 12–13. 3. NOVAK, G. Papir, karton, lepenka. V Grafični materiali. Ljubljana: Naravoslovnotehniška fakulteta, Oddelek za tekstilstvo, 2004, str. 41–158. 4. GRILJ, S. Preizkusne metode karakterizacije klasičnih in recikliranih papirjev. Seminarska naloga, Univerza v Ljubljani, Naravoslovnotehniška fakulteta, 2010, 32 str. 5. COAKLEY, T., DUFFY, N., FREIBERGER, S., FRESNER, J., HOUBEN, J., KERN, H., KRENN, C., MCCARTHY, C. in RAUPENSTRAUCH, H. Energetska učinkovitost v industriji, Priročnik za dijake. IUSES, Intelligent Energy – Europe, 2012, 78 str. 6. Ločeno zbiranje odpadkov: Papir in karton. Javni holding Ljubljana, Snaga, 2011. http:// www.jhlj.si/snaga/locevanje/papir. 7. KÖNIG, S. Lastnosti recikliranih papirjev in njihova tiskovna kakovost : doktorska disertacija. Ljubljana, 2013, 146 str. 8. Recycled fiber and deinking. Edited by L. Göttsching and H. Pakarinen. Helsinki : Fapet Oy; Atlanta : TAPPI, 2000, 649 str. 9. Star papir za novo upanje 2012, Gradivo za medije. http://ebm.si/p/starpapir2012/Gradivo_za_medije_SPZNU_2012_09112012. pdf. 10. Key statistics: European pulp and paper industry 2012. Bruselj : Confederation of European Paper Industries, 2013, 32 str. 11. Uredba o odpadkih. Uradni list Republike Slovenije, Ljubljana, 2011. http://www.uradnilist.si/1/content?id=106484. 12. Direktiva Evropskega parlamenta in Sveta 94/62/ES o embalaži in odpadni embalaži. Uradni list Evropske unije, Bruselj, 1994, str. 349–362. 13. Uredba o zelenem javnem naročanju. Uradni list Republike Slovenije, Ljubljana, 2011. http://www.uradnilist.si/1/content?id=106374. 14. GREGOR-SVETEC, D., MOŽINA, K., BLAZNIK, B., URBAS, R., VRABIČ BRODNJAK, U. in GOLOB, G. Efficient paper recycling. XI Symposium on Graphic Arts, Pardubice : Faculty of Chemical technology, 2013, str. 44– 51. 50 | Pomurska obzorja 1/ 2014/ 2 15. Direktiva 2008/98/ES Evropskega parlamenta in Sveta o odpadkih in razveljavitvi nekaterih direktiv. Uradni list Evropske unije, Strasbourg, 2008, str. 312/3–312/30. 16. VILLANUEVA, A. in EDER, P. End-of-waste criteria for waste paper: Technical proposals. IPTS, Sevilla, 2011, 101 str. 17. Paper recycling: European declaration on paper recycling 2011–2015. Bruselj : European recovered paper council, 2011, 8 str. Tehnika Jože Nemec Računalnik, ki je spremenil svet POVZETEK V preteklosti so podatke obdelovali le na velikih računalnikih v računskih centrih. Ko so se sredi osemdesetih let pojavili prvi, tedaj še okorni hišni računalniki, so jih pričeli, zlasti mladi, čedalje bolj množično uporabljati. Resničen preboj pa se je zgodil prišel šele leta 1981, ko so predstavili osebni računalnik IBM-PC. Četudi ga celo v sami tovarni IBM niso pretirano cenili, pa je njegova odprta arhitektura kmalu našla posnemovalce. V nekaj letih je ta računalnik, s svojimi nasledniki, postal standard za strojno opremo. Množica programske opreme in pocenitev strojne opreme je računalnikom utrla pot v vse pore življenja. Danes se le redki zavedamo kako hiter razvoj je ta računalnik dosegel na področju strojne in programske opreme. Zato je pogled na zgodovino razvoja osebnih računalnikov zanimiv in poučen. Ključne besede: osebni računalnik, strojna oprema, programska oprema. V sedemdesetih in osemdesetih letih preteklega stoletja so podjetja obdelovala podatke z velikimi računalniškimi sistemi, nameščenimi v računskih centrih. Prostor, ki so ga zasedali računalniki, je običajno presegal sto kvadratnih metrov. Računalnike so lahko programirali in upravljali le redki, posebej izšolani, strokovnjaki. Običajnim uporabnikom računalniških uslug pa je bil proces obdelave tuj. Integrirana vezja, ki so sprva vsebovala le nekaj tranzistorjev, so postajala vedno kompleksnejša. Leta 1971 je podjetje Intel poslal na tržišče prvi mikroprocesor 4004. Tudi drugi izdelovalci polvodniških elementov kot so Zylog, Texas Istruments, Siemens, so razvili svoje tipe mikroprocesorjev. Vsa ta vezja so kmalu uporabljali pri vodenju raznih obdelovalnih strojev. Svoje mesto pa so našla tudi v napravah, ki so zahtevale primerno računalniško upravljanje. V tem času so nastali tudi prvi osebni računalniki. Na tržišču so bili računalniki različnih oblik in sposobnosti. Aprila 1976 so Steve Jobs, Steve Wozniak in Ronald Wayne pričeli prodajati matične plošče za računalnik Apple I v kitu. Podjetje Heathkit je leta 1978 poslalo na tržišče računalnik H8. Tudi tega ste dobili v kitu, kar pomeni, da ste ga morali sami sestaviti. Oba računalnika sta stala okoli 600 $. Za »malenkosti« kot so tipkovnica, monitor in kasetofon, ki je služil kot zunanji pomnilnik, je moral kupec poskrbeti sam. Leta 1980 je podjetje Commodore Busines Machines pričelo prodajati računalnik Vic-20, ki je bil predhodnik pri nas popularnega računalnika Commodore C64. Istega leta je podjetje Sinclair Research Ltd pričelo izdelovati računalnik ZX Spectrum. Podobne računalnike so izdelovala še številna podjetja kot Atari, Oric in druga. Za vse te računalnike je bilo značilno, da so bila celotna vezja nameščena na dnu tipkovnice. Za monitor so uporabljali televizorje, za zunanje pomnilnike pa kasetofone. Zaradi njihove cenene izvedbe je bilo zanje na tržišču veliko igric, tako da so bili ti računalniki neke vrste predhodniki igralnih konzol. Ob teh računalnikih, ki so bili namenjeni predvsem domači rabi, pa je nastalo tudi veliko profesionalnih osebnih računalnikov. Med take sodijo Datamaster tovarne IBM iz leta 1979, TFC 3450 tovarne Fujitsu iz leta 1981, računalniki Rainbow, DECmate, Professional 325 in 350 tovarne Digital Equipment Corporation, ter številni računalniki HewlettPackard, podjetja, ki v tem obdobju izdelovalo predvsem elektronske merilne instrumente. Za vse te računalnike je bilo značilno, da so uporabljali različne operacijske sisteme, da so bili programi napisani za določen računalniški sistem in zato na drugih niso delovali. Ob tako pisani ponudbi na tržišču je, tedaj daleč največji proizvajalec računalnikov, IBM predstavil svoj osebni računalnik IBM Personal computer. Ta računalnik je bil poznan tudi po kratici IBM-PC. Z njim je IBM želel prodreti na trg osebnih računalnikov, hkrati pa bi naj te naprave služile tudi kot inteligentni terminali za njegove velike sisteme. Razvoj tega računalnika je bil kočan v enem letu, čeprav so pri IBM-u za razvoj takih sistemov običajno porabili okoli štiri leta. Slika 1. Osebni računalnik IBM-PC. Pri razvoju tega osebnega računalnika je IBM sprejel kar nekaj, zanj nenavadnih, odločitev. Četudi je imel razvit svoj mikroprocesor, je uporabil Intelov mikroprocesor 8088. V času, ko so bile proizvodne zmogljivosti integriranih vezij omejene, je lahko le Intel zagotovil tekočo dobavo svojih mikroprocesorjev. Pomurska obzorja 1/ 2014/ 2 | 51 Jože NEMEC: RAČUNALNIK, KI JE SPREMENIL SVET Prav tako so bila tudi vsa ostala integrirana vezja in elektronski elementi izdelani pri zunanjih dobaviteljih. Vse dotlej je skušal IBM vgrajevati v računalnike predvsem komponente, ki jih je izdeloval v svojih tovarnah. Ob predstavitvi računalnika je objavil tudi njegovo podrobno tehnično shemo, ki je ni avtorsko zaščitil. S tem je skušal spodbuditi druge proizvajalce, da bi proizvajali dodatke. Tudi za te dodatke ni zahteval licenčnine. Na področju programske opreme je IBM objavil in avtorsko zaščitil BIOS. To je rezidenčni del programske opreme, ki služi kot vmesnik med strojno opremo in operacijskim sistemom. Računalnik je lahko uporabljal različne operacijske sisteme. Najpreprostejši je bil IBM Basic. Ta je bil namenjen predvsem šolam, za zahtevnejše naloge pa ni bil uporaben. Pri profesionalni uporabi ste lahko izbirali med PC-DOS 1.0, podjetja Microsoft, CP/M-86 podjetja Digital Research in UCSD p-System razvit na University of California iz San Diega. Za PC-DOS se je velika večina uporabnikov odločila predvsem zaradi cene, saj je bil ta operacijski sistem tudi do desetkrat cenejši od ostalih dveh. Vse do prihoda računalnika IBM-PC je IBM prodajal svoje izdelke le preko lastne prodajne mreže. V primeru tega računalnika pa je IBM spremenil svojo prodajno strategijo. Tako ste lahko računalnik IBM-PC v Ameriki kupili pri vseh večjih trgovskih verigah, ki so prodajale elektronsko opremo. segala od preprostih opravil, kot so merjenje napetosti, pa vse do vmesnikov za različne kompleksne regulacijske sisteme. Istočasno je nekaj podjetij pričelo s kloniranjem BIOSa. Ko je bilo to opravilo uspešno končano, so se lahko pojavili tudi prvi PC kloni. Prvi računalnik, ki je imel praktično enake lastnosti kot IBM-PC, je izdelalo podjetje Columbia Data Products. Računalnik MPC 1600 je bil predstavljen junija 1982. leta, kar je manj kot leto dni po predstavitvi IBM-ovega originala. Temu so kmalu sledili tudi drugi proizvajalci in na tržišču se je pojavilo več novih računalnikov. Med njimi velja omeniti računalnik Compaq Portable podjetja Compaq. Ta prvi prenosni računalnik je tehtal nekaj manj kot 13 kg. Čeprav ni imel akumulatorskega napajanja, je bil začetnik bodočih avtonomnih prenosnih naprav. Preglednica 1. Osnovne lastnosti računalnika IIMB-PC. Procesor Intel 8088 Takt 4,7 MHz Pomnilnik 16 - 640 kByte Zunanji pomnilnik Grafika Operacijski sistem Najcenejša konfiguracija Kasetofon Disketa 360 kByte 80 × 25 znakov 640 × 200 točk IBM Basic PC DOS 1.0 CP/M-86 UCSD p-System 1565 US $ Slika 2. Compaq portable – prvi prenosni PC. Hitra uveljavitev računalnika IBM-PC in njegovih klonov je vzpodbudila razvijalce programske opreme, da so pričeli zanje pisati programe. To je spodbudilo opuščanje proizvodnje računalnikov z drugimi operacijskimi sistemi. Praktično edina izjema, ki se je obdržala na tržišču, so računalniki tovarne Apple. Za to podjetje je značilno, da ni odprlo, tako kot IBM, svoje arhitekture. S pametno strategijo je osvojilo določen segment tržišča in danes sodi med najuspešnejše proizvajalce osebnih računalnikov. Pocenitev diskovnih pogonov, prihod sposobnejših procesorjev in drugih elektronskih elementov, je vplivalo tudi na zgradbe računalnikov. IBM je na te spremembe reagiral z spremenjenimi modeli PC-XT in PC-AT. Osnovne lastnosti teh modelov so podane v preglednici 2. Oglejmo si sedaj osnovne lastnost računalnika IBM-PC. Le te so prikazane v preglednici 1. V osnovni konfiguraciji je imel le 16kByte centralnega pomnilnika, za zunanji pomnilnik je uporabljal kasetofon, za monitor pa televizijski Preglednica 2. Primerjava modelov IBM-ovih osebnih sprejemnik. Četudi je bila ta konfiguracija računalnikov. namenjena predvsem šolam in domači rabi, je bila njena cena 1565 dolarjev. Z upoštevanjem Model IBM PC IBM PC XT IBM PC AT inflacije in tekočega tečaja €, bi znašala ta cena Avgust 1981 Marec 1983 Avgust 1984 Predstavitev danes okoli 4500 €. Če je kupec navedeni 8088 8088 80286 Procesor konfiguraciji dodal dva disketna pogona, monitor 8/16 8/16 16 Vodilo bitov in iglični tiskalnik, se je cena povzpela na okoli 4,77 MHz 4,77 MHz 6 in 8 MHz Takt 13.000 €. Takšna konfiguracija je uspešno nadomeščala terminale na večjem računalniku. 640 × 200 640 × 350 640 × 480 Grafika Tako je IBM s svojim osebnim računalnikom max 640 KB max 640 KB max 16 MB Pomnilnik pokril vse segmente trga. 10 MByte 20 MByte Disk Odprtost računalniškega sistema IBM-PC je hitro povzročila premike na tržišču. Kmalu so se pojavile prve kartice za razširitvene reže. Njihova uporabnost je 52 | Pomurska obzorja 1/ 2014/ 2 Jože NEMEC: RAČUNALNIK, KI JE SPREMENIL SVET Iz preglednice je razvidno, da so imeli kasnejši modeli že vgrajen trdi disk. Četudi so bile kapacitete teh diskov, z današnjega vidika, majhne, so bile, v primerjavi z disketnimi pogoni, zelo velike. Procesorji 8088 so imeli 8 bitno zunanje vodilo, četudi so bili interno 16 bitni. Procesorji 80286 pa so imeli tudi 16 bitno zunanje vodilo. Zato so ti procesorji lahko uporabljali večji pomnilnik. Bili so tudi hitrejši, zato so običajno imeli tudi tipko, s katero se je zmanjšala hitrost delovanja. To zmanjšanje je bilo potrebno uporabiti takrat, ko je višja hitrost povzročila motnje pri delovanju starejših programov. Leta 1985 je podjetje Intel predstavilo procesor 80386. Ta 32 bitni procesor je imel 275.000 tranzistorjev in ves nabor ukazov procesorjev 8088 in 80286. Ko je ta procesor proti koncu leta 1986 prišel na tržišče, je bilo podjetje Compaq prvo, ki ga je uporabilo pri svojih računalnikih. Ti so lahko delali s taktom 33 MHz in imeli trde diske s kapaciteto 30 MByte in več. Imeli so disketno enoto s 5¼ colskimi disketami. V računalnik ste lahko dodali še disketno enoto s 3½ colskimi disketami in tračno enoto. Oglejmo si še razvoj operacijskih sistemov za IBM-PC. Predhodnik prvega operacijskega sistema je QDOS podjetja Seattle Computer Products, ki ga je predstavilo junija leta 1979. Njegov avtor Paterson ga je poimenoval Quick and dirty operating system. S tem je povedal, da je bil napisan v naglici in površno. V času nastanka IBM-ovega osebnega računalnika Microsoft ni imel ustreznega operacijskega sistema, zato je kupil julija 1981 QDOS in ga rahlo dodelal. IBM je ta operacijski sistem ponujal kot PC-DOS. Maja leta 1982 pa je Microsoft poslal na tržišče OEM verzijo tega operacijskega sistema z oznako MS-DOS 1.25. Kasneje so prihajale na tržišče nove različice tega operacijskega sistema. V njih so bile odpravljene določene napake, istočasno pa so omogočale tudi uporabo dodatnih komponent kot so 3½ colske diskete, trdi diski in podobno. Prva verzija operacijskega sistema Windows je bila predstavljena novembra 1985. Bila je le razširitev MS-DOS-a z novim grafičnim vmesnikom. Prva, zares uporabna verzija tega programa, Windows 3.11 pa je prišla na tržišče aprila 1991. Leto za tem je borzna vrednost Microsofta presegla borzno vrednost IBM-a. Naslednje verzije so imele za oznako letnice predstavitev kot Windows 95 in Windows 98. Zadnja verzija Windows 8 pa ponovno uporablja številčno oznako. Na vse večje izzive proizvajalcev klonov je IBM aprila 1987 odgovoril s predstavitvijo nove družine računalnikov. Zanje je uporabil naziv IBM Personal System oziroma IBM-PS/2. Ta družina je imela vrsto modelov, ki so se razlikovali predvsem po procesorjih. Osnovni model je imel še 8/16 bitni procesor 8086, medtem ko so imeli ostali modeli 16 bitni procesor 80286 in 32 bitni procesor 80386. Ostale osnovne lastnosti so prikazane v preglednici 3. Preglednica 3. Osnovne lastnosti sistemov IBM-PS/2. Pomnilnik 1 - 16 MByte Takt 8 - 90 MHz Grafika 1024 × 768 Zunanji pomnilnik Disk ≥ 40 MB Disketa 3½ col V verzijah z zmogljivejšimi procesorji je skušal IBM uveljaviti novo vodilo, ki je kasneje postavilo standard za vodilo ISA. Uvedel je še druge lastnosti kot so novi priključki za tipkovnico in miško. Poglavitna sprememba IBM-ove strategije pa je bila, da je za uporabo novega vodila zahteval licenčnino. IBM je skušal uveljaviti tudi nov operacijski sistem OS/2. Vse te novosti niso naletele na pričakovan sprejem. Namesto pojma kompatibilno z IBM-PC, se je vedno bolj uveljavljal pojem MSDOS oziroma Windows kompatibilen računalnik. Leta 2004 je IBM svoj oddelek za proizvodnjo osebnih računalnikov prodal podjetju Lenovo. To podjetje je še nekaj let proizvajalo računalnike z oznako IBM. Danes je IBM na področju informacijske tehnologije še vedno eno izmed vodilnih podjetij na svetu. To velja tako za področje programske, kot tudi strojne opreme. Pojav osebnih računalnikov pa je vplival tudi na druga podjetja, ki so se v preteklosti ukvarjala z informacijsko tehnologijo. Tako se je Hewlett-Packard, iz podjetja za proizvodnjo merilnih inštrumentov, spremenil v podjetje, ki se ukvarja izključno z računalniško tehnologijo. V preteklosti je HP tudi kupil tako zveneča podjetja kot sta Compaq in DEC. Pojavila so se številna podjetja, ki ponujajo programske rešitve za osebne računalnike. Številna podjetja razvijajo nove dodatke in naprave za osebne računalnike. Vse to pa omogoča, da so osebni računalniki postali naprave, ki jih srečamo skoraj na vsakem koraku. In četudi se počasi bližamo tehnološki meji, ki jo omogočajo polvodniki, bomo priča še številnim inovacijam in njihovem vplivu na vsakdanje življenje. Literatura 1. IBM 5150 Technical Reference 6025005, IBM, avgust 1981 2. Compaq Deskpro 386/33 Personal Computer, Features/Specifications, Compaq Computer Corporation, maj 1989 3. The IBM Personal Computer, First Impressions, Byte, oktober 1981 4. http://en.wikipedia.org/wiki/IBM_Personal_Computer 5. http://en.wikipedia.org/wiki/IBM_PC_compatible 6. http://www03.ibm.com/ibm/history/exhibits/pc25/pc25_fact.html 7. http://en.wikipedia.org/wiki/MS-DOS 8. http://sl.wikipedia.org/wiki/Microsoft_Windows Pomurska obzorja 1/ 2014/ 2 | 53 Tehnika Borut Žalik1*, Niko Lukač1, Domen Mongus1 Algoritmi nad velikimi količinami podatkov LiDAR POVZETEK V članku predstavimo postopke obdelave velikih količin nestrukturiranih geometrijskih podatkov, pridobljenih s prebirniki LiDAR. Najprej opišemo postopke klasifikacije točk LiDAR v točke, ki ležijo na Zemljinem površju. Opravimo kratek pregled obstoječih metod nato pa na kratko predstavimo našo metodo, temelječo na konceptih matematične morfologije. V nadaljevanju predstavimo zanimivo aplikacijo klasificiranih podatkov LiDAR za določitev sončevega potenciala streh v urbanih področjih. Podamo ključne parametre metode in predstavimo podatke za določitev sončevega potenciala in pridobljene električne energije iz treh različni h tipov celic za občino Beltince. Ključne besede: računalništvo, daljinsko zaznavanje, LiDAR, klasifikacija. 1. Uvod Zajemanje podatkov o površju Zemlje je od nekdaj vznemirjalo človeštvo, saj je temelj za naše bivanje v okolju. Ključen podatek pri opisovanju površja je položaj točk v ustreznem koordinatnem sistemu, med katerimi izmerimo še kote in razdalje. Tehnologije zajemanja in opisovanja točk na zemeljskem površju so se seveda spreminjale. Najzgodnejši postopki izpred nekaj 1000 let so uporabljali le količek in vrvico, kasneje so jih zamenjeli različni instrumenti. Ena izmed najpriljubljenejših tehnologij današnjega časa je tehnologija LiDAR (Light Detection And Ranging) [1]. Gre za tako imenovano aktivno tehnologijo daljinskega zaznavanja, ki temelji na merjenju razdalje od naprave do objekta z lasersko svetlobo. Ker se svetloba od objektov na površju odbija, nekaj žarkov doseže tudi tla v gosti vegetaciji. Za zajem zemeljskega površja pritrdijo merilnike LiDAR na letalo ali helikopter, naprava pa zatem odbere tudi do 200.000 točk v sekund in, odvisno od višine in hitrosti leta, doseže gosto nekaj 10 točk na kvadratni meter. Rezultat je ogromna količina nestrukturiranih podatkov (3D točk s prirejenimi parametri, kot so intenziteta odboja, številka odboja, čas zajema), iz katerih želimo povsem samodejno z računalniškimi algoritmi pridobiti uporabne podatke za nadaljnje obdelave. Učinkoviti algoritmi za obdelavo takšnih podatkov danes predstavljajo glavni raziskovalni izziv z množico aplikacij [2, 3]. Osnovni problem opišemo z lahkoto. Naj bo P množica n točk pi, 0 < i < n, pridobljenih iz prebirnika LiDAR. Vsako točko pi želimo uvrstiti v eno izmed podmnožic, ki pripadajo terenu, objektom na terenu (na primer hiše), 1Univerza v Mariboru, Fakulteta za elektrotehniko, računalništvo in informatiko, Laboratorij za geometrijsko modeliranje in algoritme multimedijev, / Smetanova 17, 2000 Maribor, Slovenija E-naslovi: borut.zalik@um.si; niko.lukac@uni-mb.si; domen.mongus@uni-mb.si *Avtor za korespondenco 54 | Pomurska obzorja 1/ 2014/ 2 vegetaciji (nizki ali visoki vegetaciji) ali šumu v podatkih. Postopku pravimo klasifikacija. V prispevku bomo predstavili rezultate lastnih algoritmov klasifikacije, osredotočili pa se bomo na klasifikacijo terena. Kot primer uporabe klasificiranih podatkov LiDAR bomo pokazali uporabo podatkov LiDAR za izračun fotovoltaičnega potenciala streh občine Beltinci. Članek je razdeljen v štiri poglavja. V poglavju 2 predstavimo algoritem za določanje točk terena, v tretjem poglavju opišemo aplikacijo, ki na temelju klasificiranih podatkov LiDAR izračuna fotovoltaični potencial streh objektov, v četretk poglavju pa povzamemo članek. 2. Klasifikacija točk terena iz oblaka točk LiDAR Klasifikacija točk terena velja za enega izmed temeljnih problemov obdelave podatkov LiDAR in je predmet številnih raziskav od samega začetka razvoja tehnologije [4]. Čeprav vse poznane metode temeljijo na dejstvo, da so točke terena nižje od točk objektov nad njim, se specifike posamičnih pristopov med seboj pomembno razlikujejo. Poznane metode, imenovane tudi metode filtriranja podatkov LiDAR, v grobem delimo na [5,6]:  napovedovalne metode filtriranja, ki se zanašajo na grobo aproksimativno površje, na osnovi katerega izvedejo linearno napovedovanje,  metode, ki temeljijo na ocenitvi naklona in izhajajo iz predpostavke, da je naklon točk objektov do okoliškega terena občutno večji, kot je naklon med samimi točkami terena ter  metode matematične morfologije, ki točke LiDAR najprej razporedijo v mrežo enako velikih celic, nad to mrežo pa izvajajo postopke računalniškega vida razvite za obdelavo slik. V nadaljevanju predstavimo. vsak pristop nekoliko podrobneje Borut ŽALIK, Niko LUKAČ, Domen MONGUS: ALGORITMI NAD VELIKIMI KOLIČINAMI… 2.1. Napovedovalne metode filtriranja podatkov LiDAR Osnovna ideja, na kateri temeljijo napovedovalne metode filtriranja podatkov LiDAR, je bila predstavljena v [7] in jo povzema slika 1. Prvi korak filtriranja z linearnim napovedovanjem predstavlja izračun aproksimativne funkcije, ki daje grobo oceno oblike terena. V ta namen so avtorji v začetku predlagali uporabo linearne aproksimacije z metodo uteženih najmanjših kvadratov. V prvi iteraciji določimo enake utežne vrednosti vsem vhodne točke LiDAR, skozi nadaljnje iteracije pa jih posodabljamo z linearno utežno funkcijo. Utežna funkcija sledi predpostavki, da točke nad aproksimativnim površjem pripadajo objektom in jim zato zmanjša utežne vrednosti, medtem ko utežne vrednost točkam pod aproksimativnim površjem poveča. Točkam, katerih višinske razlike so nad vnaprej določeno pragovno vrednostjo, določi utežno vrednost 0, s čimer jih izloči iz nadaljnjega postopka. S postopkom zaključimo, ko sistem preide v stabilno stanje in se utežne vrednosti točk ne spreminjajo več oziroma ko dosežemo vnaprej določenem številu iteracij. aproksimacija daje zelo zglajeno oceno površja, ti pristopi pogosto izločijo detajle terena (na primer ostre robove, prepade ter nasipe) in pomotoma ohranijo majhne in nizke objekte (na primer nizko vegetacija, ograje in avtomogile) [4,5,6]. 2.2. Metode filtriranja podatkov LiDAR glede na naklon Metode, ki so zmožne odpraviti nekatere pomanjkljivosti pristopov z linearnim napovedovanjem (vendar ne vseh), temeljijo na oceni naklon med točkami in njihovo okolico. Enega izmed prvih takšnih pristopov je predlagal Vosselman [10]- Izhajal je iz predpostavke, da se naklon terena ne dvigne nad določeno pragovno vrednost. Takšna predpostavka daje možnost pragovnega filtriranja točk glede na gradient, ki ga te tvorijo s svojo soseščino. V ta namen je predlagal pragovno funkcijo v obliki stožca, ki ga določa polmer osnovnega kroga v ravnini xy ter višina glede na koordinato z. Hitro lahko opazimo, da polmer kroga določa preiskovalno okolico, višina pri danem polmeru pa pragovno vrednost naklona. Postopek filtriranja nato izvedemo tako, da vrh stožca premaknemo v dane točke ter preverimo, ali katera izmed točk iz njegove okolice leži pod Slika 1. Filtriranje podatkov LiDAR na osnovi linearnega napovedovanja, kjer (a) vhodnim podatkom LiDAR najprej določimo (b) začetni približek terena ter izvedemo (c) pragovno filtriranje, s katerim izločimo točke, ki ne pripadajo terenu. Točkam nato posodobimo utežne vrednosti. Z iterativnim ponavljanjem teh dveh korakov dosežemo (d) klasifikacijo točk terena. Očitno je, da je natančnost tega pristopa v veliki meri odvisna od izbrane pragovne vrednosti in uspešnosti prilagajanja utežnih vrednosti točk. Avtorja v svoje kasnejšem delu podrobneje analizirata številne definicije utežnih funkcij [8] in pri tem pokažeta, da se optimalna pragovna vrednost spreminja glede na število iteracij. Kot rešitev predlagata vektor pragovnih vrednosti. Vseeno pa sorodna dela izpostavlja težave metode pri klasifikaciji točk terena s spremenljivim naklonom. Razgibani tereni namreč zahtevajo dinamično prilagajanje utežnih funkcij. Lee in Younan [9] zato predlagata prilagodljivo filtriranje normaliziranih višinskih razlik, kjer uvedeta številne prilagoditvene parametre, na primer parameter faktor zakasnitve in stopnja filtriranja. Čeprav tako občutno izboljšata natančnost pristopa predvsem v najzahtevnejših primerih, tudi ta metoda ni zmožna odpraviti poglavitnih pomanjkljivosti pristopov, ki temeljijo na grobi aproksimaciji površja. Ker takšna njim. Če pod stožcem najdemo eno ali več točk iz okolice, potem dano točko označimo kot točko objekta in jo izločimo. Kadar pa pod stožcem ne leži nobena točka, lahko dano točko označimo kot del terena (glej sliko 2a). Uspešnost tega pristopa je v prvi vrsti odvisna od pravilno določenih vrednosti naklona, oziroma višine stožca in njegovega polmera. Ker osnovna različica metode predvideva konstanten naklona, tudi ta metoda ni primerna za obdelavo razgibanih terenov s spremenljivim naklonom. Še zlasti kadar teren vsebuje nezvezne predele, podana definicija pragovne funkcije pogosto povzroči njihovo prefiltriranje in vodi do zglajenega digitalnega modela reliefa, kot to prikazuje slika 2. Pomurska obzorja 1/ 2014/ 2 | 55 Borut ŽALIK, Niko LUKAČ, Domen MONGUS: ALGORITMI NAD VELIKIMI KOLIČINAMI… Slika 2. Filtriranje podatkov LiDAR glede na naklon, kjer je (a) pragovna funkcija določena v obliki stožca premeščenega v točko, ki jo trenutno filtriramo. Dano točko označimo kot del terena, kadar se pod stožcem ne nahaja nobena točka iz okolice. V nasprotnem primeru točko označimo kot del objekta in jo izločimo. Takšen pristop vodi do (b) prefiltriranja nezveznih predelov terena in posledično (c) zglajenega digitalnega modela reliefa. Eden izmed prvih pristopov k odpravi te pomanjkljivosti je predstavil Sithole [11]. V ta namen je podrobneje analiziral več tehnik ocenitve lokalnega naklona terena, iz katerih lahko izpeljemo učinkovitejšo definicijo pragovne funkcije. Čeprav je tako uspel občutno izboljšati odpornost filtra na lokalne variacije terena, tega ni uspel povsem odpraviti. Zaradi simetrične definicije pragovne funkcije predlagani pristop še vedno ni zmožen uspešno filtrirati nezveznih in nesimetričnih predelov terena. Posebno težavo predstavljajo točke lokalnih ekstremov in prevojev, kjer se gradienti v različnih smereh občutno razlikujejo. Možno rešitev tega problema sta nedavno predstavila Wang in Tseng [12], ki sta predlagala tako imenovani dvosmerni prilagodljiv filter, ki temelji na delitvi pragovne funkcije (stožca) na dve polovici. Vsak izmed filtrov je tako zmožen ohraniti nezveznosti terena v nasprotni smeri, rezultat filtriranja pa tako predstavlja unija točk, ohranjenih z enim ali drugim filtrom. Čeprav takšen pristop vodi do ohranitve detajlov terena, pa občutno poslabša izločanje objektov, še zlasti kadar so ti veliki in nizki (na primer skladišča). Tovrstni pristopi prav tako ne nudijo najvišje stopnje natančnosti pri filtriranju točk nizke vegetacije na strmih pobočjih. Točke nizke vegetacije imajo namreč relativno nizek naklon do svoje okolice, medtem ko je naklon terena lahko visok. Posledično to pomeni, da začetna predpostavka, na kateri temeljijo predstavljeni postopki filtriranja, ne drži več [4,5,6]. Zato nedavno predstavljena študija Werbroucka in sodelavcev [13] predlaga uporabniško nadzorovano pragovno filtriranje naklonov z vizualno interpretacijo podatkov. Takšen pristop seveda zahteva zelo intenzivno uporabniško interakcijo, a je to edini način, s katerim lahko dosežemo visoko stopnjo natančnosti na osnovi ocenitve naklonov. 2.3. Metode matematične morfologije Čeprav je bila prva metoda filtriranja podatkov LiDAR s koncepti matematične morfologije predstavljena že leta 1996, so ti pristopi (kot bomo predstavili v nadaljevanju) komaj v zadnjem času prevzeli vodilno vlogo na tem področju tako iz vidika natančnosti kot tudi računske učinkovitost. Vsem tovrstnim postopkom je skupno to, da vhodni oblak točk najprej razporedijo v mrežo enako velikih celic, nad njo pa izvajajo operacije, izpeljane iz morfološkega širjenja in krčenja [14]. Slika 3. Primer izravnavanja detajlov na (a) prečnem prerezu vhodnih podatkov LiDAR, (b) razporejenih v mrežo enako velikih celic z uporabo (c) morfološkega krčenja ter (d) širjenja. 56 | Pomurska obzorja 1/ 2014/ 2 Borut ŽALIK, Niko LUKAČ, Domen MONGUS: ALGORITMI NAD VELIKIMI KOLIČINAMI… Slika 3 prikazuje primer uporabe morfološkega širjenja in krčenja, ki vodi do najpogosteje uporabljenega operatorja, ki ga imenujemo morfološko odpiranje. Na sliki 3a lahko vidimo prečni prerez vhodnega oblaka točk, ki ga razporedimo v mrežo enako velikih celic glede na ravnino xy. Slika 3b prikazuje prečni prerez mreže, kjer je vsaka celica predstavljena v obliki stolpca in vsebuje natanko eno točko. Operator morfološkega krčenje (prikazan na sliki 3c) določi novo višinsko vrednost točke tako, da ji pripiše najnižjo višinsko vrednost točke v njeni okolici, določeni s strukturnim elementom W. Nasprotno pa operator širjenja (prikazan na sliki 3d) določi novo višinsko vrednost točke glede na najvišjo vrednost točke iz njene okolice. Če krčenje in širjenje izvedemo zaporedno, izvedemo izravnavo izstopajočih visokih predelov, kot to prikazuje slika 3d. Tej operacijo pravimo morfološko odpiranje, velikost izravnanih detajlov pa določa strukturni element. Prvi postopek filtriranja podatkov LiDAR z uporabo morfološkega odpiranja je predlagal Kilian s sodelavci [15]. Postopek iterativno povečuje strukturni element ter tako postopno izravna vedno večje objekte. V vsakem koraku izračuna višinsko razliko med točkami pred in po uporabi morfološkega odpiranja ter zabeleži velikost filtra pri katerem je ta prerastla določeno pragovno vrednost. To vrednost uporabi za izračun utežne vrednosti aproksimativne funkcije s katero določi digitalni model reliefa. Čeprav je postopek relativno naiven, nakazuje poglavitno težavo, s katero se srečamo pri uporabi morfoloških filtrov za klasifikacijo podatkov LiDAR. Ker je velikost izravnanih detajlov odvisna od velikosti strukturnega elementa, glavni izziv teh pristopov predstavlja ustrezna definicija strukturnega elementa. Majhen strukturni element namreč ni zmožen odstraniti velikih objektov (na primer hiš). Nasprotno pa velik strukturni element poškoduje dejale terena pri odstranjevanju velikih objektov. Zardi tega, podobno kot v primeru predstavljene metode, tudi številne kasnejše metode temeljijo na iterativnem pristopu s postopnim spreminjanjem strukturnega elementa [16]. Chen s sodelavci [16] je predstavil metodo prilagojeno za terene s spremenljivim naklonom, ki temelji na predpostavki, da objekti povzročijo opazno večje spremembe naklona kot so spremembe v naklonu terena. Filtriranje točk tako dosežemo s pragovnimi vrednostmi, ki so določene glede na gradient v dani točki. Zaradi tega pa metoda ni zmožna ohraniti ostrih nezveznih predelov terena, na primer skal in nasipov[5]. Izboljšavo je predlagal Qi s sodelavci [17], ki izhaja iz predpostavke, da se višinske vrednosti robov terena spreminjajo, medtem ko imajo objekti (npr. zgradbe) enako visoke robove. Očitno je, da predstavljena metoda ni možna izločiti objekte kompleksnejših oblik. Pristop, ki se oddaljuje od analize gradientov za filtriranje podatkov LiDAR in temelji na postopnem zmanjševanju strukturnega elementa, sta predlagala Mongus in Žalik [4]. Tako imenovani hierarhični morfološki pristop temelji na postopen približevanju interpolacijske ploskve dejanskemu terenu s postopnim vključevanjem novih točk. V vsakem koraku se v ta namen izvede morfološko odpiranje višinskih razlik med točkami in približka terena ter izloči točke, ki pri tem presežejo pragovno vrednost, določeno s statistično analizo ustvarjenih višinskih razlik. Izločene točke interpoliramo in v naslednji približek terena vključimo le najnižjo točko vsaki celici podporne mreže. Predstavljen postopek je nedavno izboljšal Chen s sodelavci [19], ki je odpravil potrebo po statistični analizi za določanje pragovnih vrednsoti, Pingel s sodelavci [18] pa je predlagal postopek, ki uporablja interpolacijsko funkcijo zgolj za določanje pragovnih vrednosti pri postopnem zmanjševanju strukturnega elementa. Čeprav te metode danes veljajo za najbolj natančne, so iterativne in posledično računsko zahtevne. Mongus in Žalik [20] sta zato predlagala hierarhični pristop, ki namesto klasičnih morfoloških operatorjev uporablja tako imenovane atributne filtre. To so filtri, ki namesto vnaprej določenega strukturnega elementa filtrirajo podatke zgolj glede na določeno njegovo lastnosti. Avtorja sta v ta namen uporabila površino, postopek pa s postopnim povečevanjem površine filtrira doseže hierarhično dekompozicijo podatkov v enem samem prehodu. Čeprav je natančnost predlaganega postopka primerljiva s sorodnimi metodami, pa ta pristop ne omogoča odstranjevanja pripetih objektov, na primer mostov. 3. Določanje fotovoltaičnega potenciala streh iz klasificiranih točk LiDAR Ena izmed zelo zanimivih aplikacij uporabe podatkov LiDAR je določanje sončnega potenciala streh, ki temelji na fotovoltaičnem pojavu (ang. photovoltaic, PV) [22]. Ker tehnologija pomembno vpliva na zmanjševanje onesnaževanja ozračja in toplo-grednih emisij, je vse bolj razširjena inštalacija sistemov PV v urbanih območjih. Pri tem je eno iz med ključnih vprašanj iskanje dovolj primernih površin na strehah stavb za inštalacijo tovrstnih sistemov, da se investicija v celoti povrne. Primernost površja je predvsem odvisna od geografske lokacije, naklona, orientacije, lokalne klime ter vpliva senčenja iz okolice. Sončni potencial predstavlja najbolj primerno metriko, ki upošteva vse naštete faktorje vpliva. Izračuna se s pomočjo povprečnega dnevnega ali kumulativnega obsevanja nad dano površje skozi celotno leto. Sončno obsevanje se lahko izračuna s pomočjo direktne, difuzne in odbojne komponente sončnega sevanja. Difuzno obsevanje je posledica perturbacij direktnega obsevanja (npr. atmosferskega sipanja, oblakov, in onesnaženost zraka) med Soncem in danim površjem na Zemlji. Direktno obsevanje prav tako ne doseže površij, ki so osenčena. Praviloma prejmejo ta površja le difuzno obsevanje iz vidnega dela neba. Odbojno obsevanje se pojavi, kadar se difuzno ali direktno obsevanje odbije od nekega površja, kjer je potrebno upoštevati odbojni koeficient materiala. Vse komponente obsevanja se v praksi izračunajo kot dobri približki, saj je simulacija lokalnega podnebja oziroma dnevnega vremena za eno leto vnaprej praktično nemogoča. Dober približek je uporaba večletnih meritev globalnega in difuznega obsevanja s piranometrom. V zadnjem desetletju je prav tako možno uporabiti satelitske podatke o spektru obsevanja [23]. Čeprav so lokalne meritve na površju bolj natančne, imajo satelitske meritve prednost na hribovitimi površji, kjer bi interpolacija lokalnih meritev imela večjo napako. Podatki LiDAR omogočajo odlično analizo topografskih značilnosti površja (t.j. naklona in orientacije) ter analizo senčenja okolice. V zadnjih letih je bilo razvitih več rešitev, ki uporabljajo podatke LiDAR za izračun sončnega potenciala. Vöegtle s Pomurska obzorja 1/ 2014/ 2 | 57 Borut ŽALIK, Niko LUKAČ, Domen MONGUS: ALGORITMI NAD VELIKIMI KOLIČINAMI… sodelavci [24] je analiziral sončni potencial nad avtomatsko določenimi ravnimi deli streh različnih stavb. Kassner in sodelavci [25] so predstavili podatke LiDAR v 2.5D mreži, kjer so nato nad zamaskiranimi površji streh izračunali sončni potencial. Jochem je s sodelavci [26] razvil natančnejšo metodo za izračun sončnega potenciala, kjer je upošteval senčenje delov streh iz okolice. Za izračun sončnega potenciala so uporabili algoritem r.sun [27], ki uporablja indeks CSI (ang. clear-sky index) za aproksimacijo difuznega in direktnega sončnega obsevanja. Jochem in sodelavci [28] so nato predstavili metodo za izračun sončnega potenciala, ki upošteva transparentno senčenje vegetacije. Pri tem so si pomagali z razmerjem med prvim in zadnjim odbojem laserskega žarka. Natančnost takšnega senčenja je tako odvisna od števila posnetkov na lokaciji skozi leto. Levinson je s sodelavci [29] razvil bolj natančno metodo za ugotavljanje vpliva senčenja od vegetacije s pomočjo modelom za posamezno vrsto vegetacije. Strehe stavb niso edina atraktivna površja za namestitev sistemov PV, saj so nekatere stene stavb prav tako primerne. Stene objektov lahko zelo natančno zajamemo z mobilnim laserskim prebiranjem, kar je naredil Jochema ter njegovih sodelavcev v [30]. Investitorje za izkoriščanje sončne energije zanima količina proizvedene električne energije, kar pa sončni potencial ne omogoča. Zato je nastala metrika fotovoltaičnega potenciala, ki upošteva karakteristiko polprevodnega materiala modulov PV ter razsmernika, nameščenega na sistemu PV. Hofierka in Hanuk [31] sta algoritem r.sun aplicirala v urbanem območju. Pri tem sta izračunala potencial PV za module iz kristalnega silicija. Wiginton je s sodelavci [32] določil potencial PV na večjih regijah s tehnologijo GIS ter razpoznavanjem stavb nad zajetimi aeroposnetki. Brito in sodelavci [33] so predvideli potenciala PV za del Lizbone z orodjem SolarAnalyst. Jakubiec in Reinhart [34] pa sta izračunala potencial PV nad urbanim območjem s pomočjo simulacije Daysim. Prav tako sta upoštevala temperature streh. Pred kratkim sta Nguyen in Pearce [35] ocenila PV potencial nad podatki LiDAR za polikristalne in amorfne silicijeve celice. V tem poglavju predstavimo izboljšano metodo za izračun PV potenciala, ki temelji na predhodnem delu iz sončnega potenciala [36] ter upoštevanja karakteristike različnih modulov PV in razsmernika. V podpoglavju 3.1 bomo najprej na kratko povzeli metodo za izračun sončnega potenciala ter predstavili metodo za izračun potenciala PV. Metodo smo aplicirali nad klasificiranimi podatki LiDAR občine Beltinci, kar podrobneje predstavimo v podpoglavju 3.2. W I i (t )  Ibi (t )  1  Si (t )   Idi (t )  2  m  kjer je t časovni trenutek, Ibi direktno obsevanje, Si koeficient senčenja ter Idi količina difuznega obsevanja za celico Ci. Direktno in difuzno obsevanje pridobimo iz meritev s piranometrom čim bliže dani lokaciji. Pri tem smo upoštevali vpadni kot sončnih žarkov Φi, ki ga izračunamo s pomočjo kotov βi in γi [37]. Koeficient senčenja izračunamo za vsako celico v mreži za določen časovni korak, ter sproti spreminjamo položaj Sonca z algoritmom SOLPOS [38]. Dana celica Ca je lahko senčena tudi od neke druge celice Cb v primeru, da velja pogoj za Si = 1: 1 S i (t )   0 Cb .z  l 2 (Ca , Cb )   Ca .z sicer kjer l2 definira l2-normo oz. evklidsko razdaljo, µ pa višinski korak za horizontalno razdaljo dveh sosednjih celic na poti med Ca in položajem Sonca. Podatki LiDAR običajno ne obsegajo celotnega področja, ki lahko vpliva na senčenje streh. Zato se lahko uporabimo tudi večločljivostno senčenje za doseganje večje natančnosti. Za to smo uporabili podatke DMR širše okolice z nižjo ločljivostjo za izračun senc okolice po istem postopku kot za visokoločljivostne podatke LiDAR [36]. V metodi smo prav tako upoštevali časovnoodvisno senčenje vegetacije, ki temelji na podatkih LAI (ang. Leaf Area Index) [39]. LAI ponazarja količino oz. gostoto vegetacije projicirane na vodoravno ravnino pod vegetacijo. V primeru senčenja od vegetacije posodobimo koeficient Si: Si (t )  max  Si (t ),1  e KL  kjer je K faktor prepustnosti svetlobe, L pa razmerje med ploščino projicirane krošnje na ravnino ter ploščino celotni ravnini pod vegetacijo. Na sliki 4 vidimo difuzno in direktno obsevanje ter primer senčenja od vegetacije. 3.1. Metoda za izračun sončnega in fotovoltaičnega potenciala Klasificirane podatke LiDAR vstavimo v 2.5D mrežo Ψ iz n celic z ločljivostjo 1m2, kar vzpostavi topologijo za ugotavljanje topografskih značilnosti, kot sta naklon βi in orientacija γi za dano celico Ci ϵ Ψ; i={1,2,…,n}. Za vsako celico v mreži izračunamo sončni potencial tako, da integriramo izračunana instantna sončna obsevanja z določenim časovnim korakom. Instantno prostorsko in časovno odvisno sončno obsevanje za dano celico je definirano kot: Slika 4. Določanje instantnega sončnega obsevanja. Vpadni kot Φi definiramo z normalo Ci.n za dano celico. Na sliki vidimo tudi vpliv senčenja od vegetacije. Metoda ne upošteva obsevanja zaradi odbojev, saj teh informacij iz podatkov LiDAR ne moremo pridobiti. Sončni 58 | Pomurska obzorja 1/ 2014/ 2 Borut ŽALIK, Niko LUKAČ, Domen MONGUS: ALGORITMI NAD VELIKIMI KOLIČINAMI… potencial nato izračunamo kot kumulativno obsevanje za neko časovno obdobje:  Wh  SPi (t )   I i (t ) dt  2  . m  Fotovoltaični potencial smo določili s pomočjo nominalne karakteristike učinkovitosti (t.j. Watt-peak [Wp]) za določen polprevodniški material (glej tabelo 1). Pri tem smo uporabili instantno sončno obsevanje, ki smo ga pomnožili s faktorjem učinkovitosti η ϵ [0, 1]. 3.2. Analiza nad občino Beltinci Metodo smo uspešno aplicirali na podatke LiDAR za občino Beltinci, kjer smo upoštevali 10-letne meritve globalnega in difuznega obsevanja (glej sliko 5) vremenske postaje v Rakičanu za časovno obdobje enega leta ter časovnim korakom 1 ure. Direktno obsevanje se je preprosto izračunalo z odštevanjem difuznega od globalnega obsevanja. Tabela 1. Faktorji učinkovitosti za tri iz med najbolj uporabljenih polprevodniških materialov. Polprevodniški material Faktor učinkovitosti (η) Amorfni silicij (A-Si) 0.16 Monokristalni silicij (M-Si) 0.11 Polikristalni silicij (P-Si) 0.08 Slika 5. Graf povprečnih meritev globalnega in difuznega obsevanja za Rakičan za zadnjih 10 let. Na sliki 6 so prikazani podatki LiDAR občine Beltinci, povprečno senčenje skozi leto, ter direktno in difuzno obsevanje. Vpliv direktnega obsevanja je posebej izrazit na površjih usmerjenih proti jugu. (a) (b) (c) (d) Slika 6. Vizualizacija centra Beltincev, kjer so (a) klasificirani podatki LiDAR, (b) povprečno senčenje skozi leto, (c) direktno obsevanje in (d) difuzno obsevanje. Pomurska obzorja 1/ 2014/ 2 | 59 Borut ŽALIK, Niko LUKAČ, Domen MONGUS: ALGORITMI NAD VELIKIMI KOLIČINAMI… Nato smo izračunali PV potencial nad danimi podatki (glej sliko 7) za tri najbolj uporabljene module PV (glej tabelo 1). Kot pričakovano, obstaja linearna odvisnost med nominalno karakteristiko učinkovitostjo modula PV ter izračunanim sončnim potencialom. Zato daje tudi monokristalni silicij (M-Si) najboljše rezultate. 4. Zaključek LiDAR je laserski merilnik, ki, če je pritrjen na satelit, letalo ali helikopter, z veliko hitrostjo, natančnostjo in gostoto zajema točke s površja terena. Medtem ko so glavni tehniški problemi samega merilnika LiDAR že uspešno rešeni, ostaja velik izziv pri avtomatski klasifikaciji zajetih točk. Točke tipično razvrstimo v točke terena, točke na objektih (na primer hiše, mostovi, ceste, železnice, daljnovodi), točke nizke, srednje in visoke vegetacije, ter šumne točke. Pričetek vsake klasifikacije je določitev točke samega terena, ki jih lahko obravnavamo kot nekakšno ozadje diskretizirane 3D slike. Klasifikaciji terena posvetimo drugo poglavje, v katerem najprej na kratko opišemo obstoječe pristope, nekoliko več časa pa se posvetimo pristopu, ki smo ga razvili v našem laboratoriju. Glede na rezultate (a) objavljene v [4 in 21] je naša metoda med najuspešnejšimi obstoječimi metodami. V nadaljevanju pokažemo eno izmed zelo zanimivih aplikacij nad klasificiranimi podatki LiDAR, to je določanje fotovoltaičnega potenciala streh na večjem urbanem področju. Predlagana metoda je trenutno najpopolnejša metoda za to nalogo, saj uporablja natančne topografske podatke samih streh (naklon, orientacija) in njihove okolice, ki jih je mogoče zajeti samo s tehnologijo LiDAR. Metoda vključuje direktno, difuzno in odbojno komponento sončnega sevanja, senčenje okoliških objektov, senčenje vegetacije glede na letni čas in večločljivostno senčenje okoliškega terena. Metoda dosega več kot 95% ujemanje glede na meritve, kot smo pokazali v [36]. Neposredna nadgradnja metode je upoštevanje karakteristik sončnih celic, s katerim lahko ocenimo količino letno proizvedene električne energije ter ocenimo čas povračila investicija. Samo metodo smo preizkusili na podatkih celotne občine Beltinci. Izdelali smo tudi spletno in mobilno inačico programske rešitve, ki na enostaven in občanom razumljiv način poda informacijo o primernosti njihove strehe za izkoriščanje sončne energije [40]. Sama aplikacija in njen pomen v lokalnem okolju bo podrobneje v svojem prispevku predstavil dr. Matej Gomboši, župan občine Beltinci. (b) (c) Slika 7. Vizualizacija centra Beltincev, kjer vidimo količino proizvedene električne energije za 1m2 za celice (a) A-Si, (b) P-Si in (c) M-Si. 60 | Pomurska obzorja 1/ 2014/ 2 Borut ŽALIK, Niko LUKAČ, Domen MONGUS: ALGORITMI NAD VELIKIMI KOLIČINAMI… International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences 2010, 38 (7B), 628–632. Zahvala Avtorji se zahvaljujejo Ministrstvu za izobrazevanje, znanost in šport ter Evropskemu skladu za regionalni razvoj (projekt DoFPoLO), Občini Beltinci za podatke LiDAR, Agenciji Republike Slovenije za okolje za podatke obsevanja ter Javni agenciji za raziskovalno dejavnost Republike Slovenije (J2-3650 in J2-5479). 13. Werbrouck, I.; Antrop, M.; Van Eetvelde, V.; Stal, C.; De Maeyer, P.; Bats, M.; Bourgeois, J.; Court-Picon, M.; Crombe, P.; De Reu, J.; De Smedt, P.; Finke, P. A.; Van Meirvenne, M.; Verniers, J.; Zwertvaegher, A. Digital elevation model generation for historical landscape analysis based on lidar data, a case study in Flanders (Belgium). Expert Systems with Applications 2011, 38 (7):8178–8185, 2011. 14. Serra, J. Image analysis and mathematical morphology; Academic Press, London, UK, 1982. 15. Kilian, J.; Haala, N.; Englich, M. Capture and evaluation of airborne laser scanner data. International Archives of Photogrammetry, Remote Sensing and Spatial Information Sciences 1996, 31 (B3), 383–388. 16. Zhang, K.; Whitman, D. Comparison of three algorithms for filtering airborne LiDAR data. Photogrammetric Engineering & Remote Sensing 2005, 71 (3), 313–324. 17. Chen, Q.; Gong, P.; Baldocchi, D.; Xie, G. Filtering airborne laser scanning data with morphological methods. Photogrammetric Engineering & Remote Sensing 2007, 73 (2), 175–185. 18. Qi, C.; Gong, P.; Baldocchi, D.; Xie, G. Filtering airborne laser scanning data with morphological methods. Photogrammetric Engineering & Remote Sensing 2007, 73 (2), 175–185. 19. Liu, X. Airborne LiDAR for DEM generation: some critical issues. Progress in Physical Geography 2008, 32 (1), 31–49. Pingel, T.J.; Clarke, K.C.; McBride, W.A. An improved simple morphological filter for the terrain classification of airborne LIDAR data. ISPRS Journal of Photogrammetry and Remote Sensing 2013, 77(1), 21–30. 20. Kraus, N.; Pfeifer, N. Determination of terrain models in wooded areas with airborne laser scanner data. ISPRS Journal of Photogrammetry and Remote Sensing 1998, 53 (4), 193–203. Chen, C.; Li, Y.; Li, W.; Dai, H. A multiresolution hierarchical classification algorithm for filtering airborne LiDAR data. ISPRS Journal of Photogrammetry and Remote Sensing 2013, 82(1), 1–9. 21. Mongus, D.; Žalik, B. Computationally efficient method for the generation of a digital terrain model from airborne LiDAR data using connected operators. IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing 2013, In press, 1–12. 22. Robinson, D.; Stone A. Solar radiation modelling in the urban context. Solar Energy 2004, 77(3), 295–309. 23. Tapiador, F. J. Assessment of renewable energy potential through satellite data and numerical models. Energy & Environmental Science 2009, 11, 1142–1161. 24. Vöegtle, T.; Steinle, E.; Tóvári, D. Airborne Laser scanning data for determination of suitable areas for photovoltaics. International Archives of Photogrammetry and Remote Sensing 2005, 36, 215–220. 25. Kassner, R.; Koppe, W.; Schüttenberg, T.; Bareth, G. Analysis of the solar potential of roofs by using official Literatura 1. 2. 3. Petrie, G.; Toth, C. K. 2008. Topographic Laser Ranging and Scanning: Principles and Processing. CRC Press: Boca Raton, USA, 2008; pp. 29–86. Brandtberg, T.; Warner, T.; Landenberger, R.; McGraw, J. Detection and analysis of individual leaf-off tree crowns in small footprint, high sampling density lidar data from the eastern deciduous forest in North America. Remote Sensing of Environment 2003, 85(3), 290-303. Brenner, C. 2000. Towards fully automatic generation of city models. The International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences 2000, 33(3), 85-92. 4. Mongus, D.; Žalik, B. Parameter-free ground filtering of LiDAR data for automatic DTM generation. ISPRS Journal of Photogrammetry and Remote Sensing 2012, 66 (1), 1–12. 5. Sithole, G.; Vosselman, G. Experimental comparison of filter algorithms for bare earth extraction from airborne laser scanning point clouds. ISPRS Journal of Photogrammetry and Remote Sensing 2004, 59 (1-2), 85– 101. 6. 7. 8. Pfeifer, N.; Reiter, T.; Briese, C.; Rieger, W. Interpolation of high quality ground models from laser scanner data in forested areas. International Archives of Photogrammetry and Remote Sensing 1999, 32 (3/W14), 31–36. 9. Lee, H. S.; Younan, N. DTM extraction of LiDAR returns via adaptive processing. IEEE Transactions on Geoscience and Remote Sensing 2003, 41 (9), 2063–2069. 10. Vosselman, G. Slope based filtering of laser altimetry data. International Archives of Photogrammetry and Remote Sensing 2000, 33, 935–942. 11. Sithole, G. Filtering of laser altimetry data using a slope adaptive filter. International Archives of Photogrammetry and Remote Sensing 2001, 34(3/W4), 203–210. 12. Wang, C.; Tseng, Y. DEM generation from airborne LiDAR data by an adaptive dualdirectional slope filter. Pomurska obzorja 1/ 2014/ 2 | 61 Borut ŽALIK, Niko LUKAČ, Domen MONGUS: ALGORITMI NAD VELIKIMI KOLIČINAMI… LiDAR data. International Archives of Photogrammetry and Remote Sensing 2008, 37, 399–403 26. Jochem, A.; Höfle, B.; Rutzinger, M.; Pfeifer, N. Automatic roof plane detection and analysis in airborne LiDAR. Sensors 2009, 9 (7), 5241–5262. 27. Šúri, M.; Hofierka, J. A new GIS-based solar radiation model and its application to photovoltaic assessments. Transactions in GIS 2004, 8 (2), 175–90. 28. Jochem, A.; Höfle, B.; Hollaus, M.; Rutzinger, M. Object detection in airborne LiDAR data for improved solar radiation modeling in urban areas. International Archives of Photogrammetry and Remote Sensing 2009, 38, 1–6. 29. Levinson, R.; Akbari, H.; Pomerantz, M.; Gupta, S. Solar access of residential rooftops in four California cities. Solar Energy 2009, 83 (12), 2120–2135. 30. Jochem, A.; Höfle, B.; Rutzinger, M. Extraction of vertical walls from mobile laser scanning data for solar potential assessment. Remote Sensing 2011, 3 (4), 650– 667. 31. Hofierka, J.; Kanuk, J. Assessment of photovoltaic potential in urban areas using open-source solar radiation tools. Renewable Energy 2009, 34 (10), 2206–2214. 32. Wiginton, L.K.; Nguyen, H.T.; Pearce, J. Quantifying rooftop solar photovoltaic potential for regional renewable energy policy. Computers, Environment and Urban Systems 2010, 34 (4), 345–357. 33. Brito, M.C.; Gomes, N.; Santos, T.; Tenedório, J.A. Photovoltaic potential in a Lisbon suburb using LiDAR data. Solar Energy 2012, 86 (1), 283–288. 34. Jakubiec, J.A.; Reinhart, C.F. A method for predicting city-wide electricity gains from photovoltaic panels based on LiDAR and GIS data combined with hourly Daysim simulations. Solar Energy 2013, 93, 127–143. 35. Nguyen, H.T.; Pearce, J.M. Automated quantification of solar photovoltaic potential in cities. International Review for Spatial Planning and Sustainable Development 2013, 1 (1), 57–70. 36. Lukač, N.; Žlaus, D.; Seme, S.; Žalik, B.; Štumberger, G. Rating of roofs’ surfaces regarding their solar potential and suitability for PV systems, based on LiDAR data. Applied Energy 2013, 102, 803–812. 37. Duffie, J.A.; Beckman, W.A. Solar Engineering of Thermal Processes. Wiley-Interscience: New Jersey, USA, 2006. 38. Reda, I.; Afshin, A. Solar position algorithm for solar radiation applications. Solar Energy 2004, 76 (5), 577– 589. 39. Yuan, H.; Dai, Y.; Xiao, Z.; Ji, D.; Shangguan, W. Reprocessing the MODIS Leaf Area Index products for land surface and climate modelling. Remote Sensing of Environment 2011, 115 (5), 1171–1187. 62 | Pomurska obzorja 1/ 2014/ 2 40. Brumen, M.; Lukač, N.; Smogavec, G.; Repnik, B.; Gomboši, M.; Žalik, B. Geografski IS za prikaz sončnega potenciala nad georeferenciranimi podatki. V: Heričko, M. (ur.); Kous, K. (ur.). V: Sodobne tehnologije in storitve: OTS 2013: zbornik osemnajste konference, Maribor, 18. in 19. junij 2013. Maribor: Fakulteta za elektrotehniko, računalništvo in informatiko, Inštitut za informatiko, 2013, str. 163-174. Tehnika Tanja Bagar* Alternativna goriva POVZETEK Kot civilizacija imamo čedalje večje potrebe po energiji, hkrati pa smo izjemno potrošniška družba, ki čedalje več kupuje in proizvaja čedalje več odpadkov. Tako iz vidika potreb in proizvodnje energije kot iz okolijskega vidika se je smiselno usmerit v proizvodnjo alternativnih goriv iz odpadnih snovi. Če pomislimo, da povprečen Slovenec proizvede 500 kg odpadkov na leto in da se od tega reciklira okrog 40%, se pojavi vprašanje kaj z ostalimi 60%. Večina teh odpadkov je seveda pristala v odlagališčih, ki kljub dobri zasnovi in vsem okoljevarstvenim ukrepom niso optimalna rešitev. Analize teh odpadkov, ki jih ne moremo snovno ponovno uporabiti so pokazale, da ima lahka frakcija zelo visoko kurilno vrednost, kar kaže na to, da bi bil to lahko izjemen energent. Ob upoštevanju okolijskih predpisov o sežigu in emisijah, predstavlja predelava neobnovljivih odpadkov v kvalitetno standardizirano alternativno gorivo (SRF) odlično rešitev, tako za okolje kot tudi za energetsko preskrbo. Če se osredotočimo na Pomursko regijo, kjer je bilo še leta 2006 odloženih 100% odpadkov, bi ponovna snovna uporaba (recikliranje in kompostiranje) in predelava v alternativna goriva, pomenila kvantni preskok v ravnanju z odpadki, saj bi odložili le 20%, 80% bi pa snovno ali energetsko uporabili. Ključne besede: alternativna goriva, termična izraba, odpadki, sežig. 1. Uvod Energetska izraba odpadkov je v svetu že zelo razvita in dodelana. Vrsta tehnologij je bila razvita, ki omogoča generacijo toplote in elektrike iz odpadnih materialov. Pri nas temu še ni tako, smo ena redkih držav, ki ne izkorišča energetske potenciala odpadkov. Zakonodaja, ki ureja področje odpadkov in predelave odpadkov v trdna goriva je zelo »birokratska« in tako je bilo v Slovenji do danes podeljeno le eno okoljsko dovoljenje za sosežig odpadkov. Vsekakor pa se celotno področje razvija in tudi mi kmalu ne bomo mogli več stati ob brzicah razvoja. Trenutna situacija je takšna, da v Sloveniji nastaja gorivo iz odpadkov ali RDF (refuse derived fuel) kot ga v tujini imenujejo, vendar ga v pri nas ne smemo oziroma nimamo kje izkoriščati, zato ga izvažamo v tujino in plačujemo za to, da ga zunaj naših mej predelajo v elektriko in toploto. Osnovna navodila ravnanja z odpadki so že dolgo znana in sprejeta (Slika 1), vendar je vprašanje koliko se v praksi izvajajo. *Zeleni Če Slovenijo primerjamo z drugimi članicami Evropske unije, daleč najbolj zaostajamo v četrtem segmentu, to je energetski izrabi odpadkov. Slika 1. Hierhija ravnanja z odpadki. rudnik Pomurja, CEROP d.o.o. E-naslov: tanja.bagar@cerop.si Tel.: +386 2 5268324; Fax: +386 2 545 93 18 Pomurska obzorja 1/ 2014/ 2 | 63 Tanja BAGAR: ALTERNATIVNA GORIVA Slika 2. Prikaz količine odpadkov na prebivalca po državah v letu 2007, porazdelitev glede na snovni krog odpadkov, vir: Eurostat, 2007. 2. Gorivo iz odpadkov Pri proizvodnji trdnega goriva iz odpadkov je potrebno upoštevati Uredbo o predelavi nenevarnih odpadkov v trdno gorivo, Uredbo o sežigu odpadkov in Uredbo o emisiji snovi v zrak iz malih in srednjih ter velikih naprav. Te uredbe predvidevajo pet razredov trdnih goriv iz nenevarnih odpadkov (Tabela 1). Ključni parametri so kurilna vrednost, klor, živo srebro, kadmij in žveplo. Na podlagi meritev se gorivo razporedi v ustrezni razred in njegova uporaba in tržna vrednost je odvisna od te uvrstitve v razred. Proizvajalci goriv se seveda trudimo, da bi alternativno gorivo iz odpadkov dosegala čim višji razred in da bi gorivo visoke kakovosti ne bi bilo več klasificirano kot odpadek, ampak kot produkt-energent. Ker lastnosti trdnega goriva, ki ga proizvajamo v podjetju CEROP d.o.o redno kontroliramo, lahko trdimo, da ima gorivo konstantno visoko kurilno vrednost (okrog 20.000 J/g) in je iz energetskega vidika primerljivo s premogom. Tabela 1. Klasifikacijski seznam trdnega goriva za razvrščanje v razrede. 64 | Pomurska obzorja 1/ 2014/ 2 Tanja BAGAR: ALTERNATIVNA GORIVA rešitev, razvijamo v podjetju produkt, trdno gorivo, ki bo zagotavljalo visoke kurilne vrednosti, hkrati pa ob ustreznem čiščenju dimnih plinov ne bo presegalo predpisanih emisijskih vrednosti. Na tem področju je v nastajanju nova uredba, ki bo predpisovala predelavo odpadkov v trdno gorivo in njegovo uporabo, vendar v tem trenutku še ni sprejeta. Ker verjamemo, da je termična izraba določenega dela odpadkov prava pot in dolgoročno najustreznejša rešitev, razvijamo takšno gorivo, katero bo možno uporabiti v obstoječih velikih kurilnih napravah. Ob ustrezni predelavi ali dodelavi že obstoječih velikih kurilnih naprav, bi lahko energetsko izrabili velik delež odpadkov. S tem bi lahko hkrati poskrbeli tako za lastno energetsko preskrbo kot tudi za najbolj smotrno odstranjevanje odpadkov na okolju prijazen način. Slika 3. Kurilne vrednosti različnih energentov in primerjava z alternativnim trdnim gorivom iz odpadkov. 3. Tehnologije termične izrabe odpadkov Ko pridelamo gorivo, ki je primerno za termično obdelavo se seveda pojavi vprašanje kam z njimi? Obstaja več tehnologij za energetsko izrabo goriva iz odpadkov, najpogostejše tehnologije, ki stojijo za procesi termične izrabe odpadkov so:  Uplinjanje (z ogljikov bogate snovi uplinimo do CO, H2 in CO2)  Termične depolimerizacija (kompleksne organske snovi se pod vplivom temperature in pritiska pretvorijo v tekočino, ki je energetsko bogata, podoben proces se dogaja med tektonskimi ploščami) 4. Zaključki Zaloge primarnih energetskih virov se trajno zmanjšujejo, potrebe po energiji pa povečujejo. Zato so alternativna goriva čedalje bolj iskana in dolgoročno tudi najbolj smiselna. Energetska izraba odpadkov je v svetu pomemben način snovne predelave odpadkov, hkrati pa nudi pomemben vir energije. Naša vizija je, da bo tudi v Sloveniji alternativno gorivo postalo pomemben in iskan energent, saj bi s tem bistveno prispevali k izpolnjevanju smernic za ravnanje z odpadki, zmanjšali izčrpavanje primarnih virov in omogočili energetsko oskrbo.  Piroliza (proces zgorevanja brez prisotnosti kisika) Literatura  Plazemsko uplinjanje (s plazemskim plamenom se spremeni organski materiale v sintetični plin in elektriko) 1. DIREKTIVA 2008/98/ES Evropskega parlamenta Evropskega sveta (19. 11. 2008) – Priloga III. 2. EN 14899:2005; Karakterizacija odpadkov-vzorčenje odpadkov – Okvirno navodilo za pripravo in uporabo načrta vzorčenja. 3. SIST EN 15442:2011; Trdno alternativno gorivo – Metode vzorčenja. 4. SIST EN 15443:2011; Trdno alternativno gorivo – Metode za pripravo laboratorijskega vzorca. 5. SIST TS CEN/TS 15413:2007; Trdno alternativno gorivo – Metode za pripravo testnega vzorca iz laboratorijskega vzorca. 6. N. Samec. Proizvodnja industrijskega goriva in okolju prijazna energija. Študija izvedljivosti projekta, 2005, pp 2-32. 7. Bosmansa A., Vanderreydtb I. , Geysenc D., Helsena L. The crucial role of Waste-to-Energy technologies in enhanced landfill mining: a technology review. Journal of Cleaner Production 55, 10-23. 8. V. S. Rotter, T. Kost, J. Winkler, B. Bilitewski, Material flow analysis of RDF Production processes, Waste management, (2004), 1005–1021. 9. N. B. Chang, Y. H. Chang, W. C. Chen, Evaluation of heat value and its prediction for refuse-derived fuel, The Science of the Total Environment, (1997), 139–148. Pri termični izrabi odpadkov je največji pomislek onesnaženje zraka. V namen ugotavljana emisij pri termični izrabi so več let spremljali emisije desetih enot a termično izrabo v Evropi in jih primerjali z evropskimi standardi. Ugotovitve so predstavljene v tabeli 2. Ta prikazuje večletno povprečje izmerjenih emisijskih vrednosti za 10 enot za termično izrabo odpadkov. Spremljali so osem emisijskih parametrov, od katerih niti en parameter ni presegel v Evropski uniji predpisanih vrednosti. Tabela 2. Klasifikacijski seznam trdnega goriva za razvrščanje v razrede. Emisije mg/Nm3 Prašni delci Povprečje 10 % od EU EU predpisi enot za sežig standarda 3,06 10 31 SO2 12,2 50 24 NOx 123 200 62 HCl 7,88 10 79 CO 26,3 50 53 Živo srebro 0,01 0,05 20 TOC 0,92 10 9 Dioksini (ng) 0,02 0,1 20 V Sloveniji zaenkrat še nimamo možnosti termične izrabe odpadkov, zato goriva, ki jih pridelamo v podjetju CEROP d.o.o izvažamo v Avstrijo. Ker to iz mnogo vidikov ni najboljša in 10. C. Montejo, C. Costa, P. Ramos, M. C. Marquez, Analysis and comparison of municipal solid waste and reject fraction as fuels for inceneration plants, Applied Thermal Engineering, (2011), 2135–2140. 11. G. Genon, E. Brizio, Perspectives and limits for cement kilns as a destination for RDF, Waste management, (2008), 2375– 2385. Pomurska obzorja 1/ 2014/ 2 | 65 Tehnika Rafael Mihalič* Simuliranje elektroenergetskega sistema – zakaj in kako POVZETEK Digitalni simulator RTDS (angl. Real Time Digital Simulator) podjetja RTDS Technologies je eden od simulatorjev, ki omogočajo izračun elektromagnetnih pojavov v realnem času, pri čemer dosega izjemno nizke korake digitalne integracije (večinoma 50 μs, za nekatere elemente močnostne elektronike, kot so na primer tiristorski pretvorniki, tudi do 1.5 μs). Posebna strojna oprema omogoča uvoz in izvoz signalov do zunanje opreme, s čimer omogoča preizkušanje dejanskih naprav v zaprti zanki z računalniškim modelom EES. Tako ima uporabnik na voljo analizo delovanja naprave same in tudi analizo vpliva delovanja te naprave na razmere v preostalem sistemu. Testirana naprava pravzaprav "ne ve" ali je priklopljena na simulacijski model ali v resničen EES. Napr avo so prvotno razvili za testiranje krmilnih modulov konverterskih postaj visokonapetostnih enosmernih prenosov. Uporaba ustrezne metode modeliranja homogenih vodov in ekstremno hitra komunikacija omogočata uporabo tehnik vzporednega procesiranja signalov, s katerim lahko pojave z ustrezno opremo izračunamo hitreje, kot bi se zgodili. Glavna težava pri razvoju takšne opreme je v doseganju dovolj velike hitrosti izračunov ob dovolj majhnem integracijskem koraku ΔT, saj manjši ΔT omogoča natančnejši izračun in zajemanje višjefrekvenčnih pojavov, po drugi strani pa zajetno povečuje potrebno število matematičnih operacij in posledično čas izračuna.. Ključne besede: elektroenergetski sistem, simulacija, simulacija v realnem času, RTDS. 1. Uvod Panta rhei (Vse teče). Ta moto, s katerim je pred približno dvema in pol tisočletjema grški filozof Heraklit (540–480 pr. n. št.) povzel svoj filozofski nazor, je osnovno vodilo pri dojemanju dinamičnih pojavov v elektrotehniki. Svet, v katerem živimo, se neprestano spreminja, nič ni stalnega; tako kot nikoli ne moremo dvakrat stopiti v isto reko, saj se ta neprestano spreminja. Tudi dogajanje v elektroenergetskem sistemu (EES) je v svojem bistvu podobno toku reke. Razmere se neprestano spreminjajo, nič ni statičnega. Stacionarna stanja, s katerimi največkrat želimo ponazoriti neko dogajanje v EES, so le bolj ali manj natančen približek dejanskih razmer, v realnosti pa prej izjema kot pravilo. V smislu inženirske prakse, katere vodilo je opis pojavov do natančnosti, ki omogoča rešitev nekega problema, je to velikokrat zadovoljivo. Vendar pa praviloma najvišje obremenitve v EES definirajo dinamični pojavi. Za upravljavce in načrtovalce EES in/ali električnih naprav je poznavanje in razumevanje dinamičnih pojavov v EES nujno. Za obravnavo dinamičnih (prehodnih) pojavov in stabilnost EES obstajajo različni matematični pristopi in orodja za izvedbo izračunov. Dinamične pojave matematično opišemo s sistemi diferencialnih enačb. V primeru, ko so spremenljivke med seboj *Fakulteta za elektrotehniko / Tržaška 25, Ljubljana E-naslov: rafael.mihalic@fe.uni-lj.si Tel.: +386 1 4768 438 66 | Pomurska obzorja 1/ 2014/ 2 linearno odvisne, ostali parametri pa konstantni (oz. lahko sistem z zadovoljivo natančnostjo na ta način predstavimo), imamo opravka z linearnimi sistemi. Matematične metode za linearne sisteme so razvite in izračun časovnih potekov veličin je bolj ali manj rutinska zadeva. Največkrat je v takih primerih moč do uporabnih rezultatov priti tudi brez simulacij in te služijo le za verifikacijo rezultatov. Vendar pa je EES nelinearen in možnost njegove linearne predstavitve je prej izjema kot pravilo. Kakor hitro je sistemska spremenljivka produkt ostalih spremenljivk (npr. moč), ne moremo več govoriti o linearnosti. Razen tega je v EES vrsta parametrov, ki so nelinearno odvisni od vrednosti sistemskih spremenljivk. Med tipične lahko štejemo npr. magnetilne krivulje transformatorjev, nasičenje generatorjev, fiksne in variabilne omejitve turbin in generatorjev, delovanje regulatorjev v najširšem smislu (turbinski, napetostni, spreminjanje odcepov transformatorjev, regulatorji elektronskih naprav ....), ki v splošnem tudi niso linearni, nadalje delovanje zaščitnih naprav, odvisnost bremen od sistemskih spremenljivk, frekvenčna odvisnost impedanc elementov EES itd. Rešitev takih nelinearnih sistemov diferencialnih enačb je mogoča samo na numeričen način, in sicer s pomočjo ustreznih računalniških programov, torej s simulacijo. Rafael MIHALIČ: SIMULIRANJE ELEKTROMAGNETSKEGA SISTEMA – ZAKAJ IN KAKO 2. Simulacija elektroenergetskih sistemov v realnem času Načeloma je pristop simulacije mogoče uporabiti pri vseh vrstah in frekvenčnih območjih prehodnih pojavov, s tem da je to ponekod edini uporaben način analize, v drugih primerih pa zopet služi za verifikacijo drugače pridobljenih rezultatov. Zgodovina izvajanja simulacij je skoraj tako dolga, kot je zgodovina EES. Ko še ni bilo mogoče izvajati numeričnih postopkov simulacije, so v ta namen uporabljali pomanjšane fizične modele EES in do rezultatov prihajali z merjenjem spremenljivk takih modelov. Obstajali so tudi hibridni modeli, nekakšen križanec med pomanjšanim modelom omrežja in računalniško tehnologijo. Ker ima v splošnem pomanjšani model večje dušenje kakor realen sistem, so s hitrimi vezji simulirali npr. »negativne upornosti« in na ta način zmanjšali dušenje modela. To tehnologijo so uporabljali zlasti za testiranje dejanske opreme v pogojih t. i. simulacije v realnem času (npr. Siemensov Schwingungsmodel). Dandanes simulacije izvajamo skoraj izključno z digitalnimi računalniki, zato pogosto srečamo tudi izraz digitalna simulacija EES. V veliki večini primerov izvajamo simulacije "off-line", kar pomeni, da je čas, kot parameter v simulaciji, neodvisen od dejanskega časa. Gre za "klasično" simuliranje z računalnikom, kjer je hitrost izračuna odvisna zgolj od hitrosti delovanja računalnika. Pri natančnem ponazarjanju ali testiranju opreme pa pristop ni tako preprost. Načeloma lahko obnašanje naprav (npr. zaščitnega releja) modeliramo z neko regulacijsko logiko. Tak pristop je povsem korekten, če imamo dovolj podatkov o napravi. Pri pojavih, ki temeljijo na nepopačenih ali malo popačenih sinusnih veličinah, v splošnem naj ne bi bilo težav. Te pa nastanejo, ko zaradi elektromagnetnih prehodnih pojavov med okvarami naprava zajema močno popačene signale, ki s sinusnimi funkcijami osnovne frekvence nimajo več veliko skupnega. Težava ni v tem, da ne bi znali modelirati logike elektronskih naprav in regulacije, temveč da je programska oprema naprav (logika) poslovna skrivnost in do podatkov praktično ni mogoče priti. Torej tega, kako se bo npr. rele obnašal pri različnih dogodkih v EES, ni mogoče v popolnosti ugotoviti, razen s t. i. simulacijo EES v realnem času. S to je namreč v realnem času mogoče simulirati dogodke v EES in signale "peljati" na napravo, izhod iz te pa zopet zajeti kot sistemsko spremenljivko v simulaciji. V realnem sistemu seveda takih preizkusov ne bi mogli narediti (npr. kratkostične okvare na 400 kV, testiranje regulatorjev naprav močnostne elektronike, zlasti HVDC …). Pri tem testirana naprava ne "ve", ali je priključena na simulator ali v resnični EES. Slika 1. Princip "real-time" simulacije v zaprti zanki Osnovno idejo simulacije v realnem času v zaprti zanki ilustrira slika 1. Simulator simulira razmere v omrežju, veličino, ki v realnem sistemu služi kot vhod v testiranec je potem potrebno prilagoditi nivoju testiranca, odziv (izhod) testiranca pa zopet uvedemo v računalnik (oz. prilagodimo nivo), kot spremenljivko, ki jo simulator upošteva v izračunih. Ključ te tehnologije predstavlja simuliranje sistema v realnem času, kar je težko uresničljiva naloga. Začetne simulacije EES v realnem času so se izvajale na analognih TNA-jih (Transient Network Analyzer). To so bili pomanjšani modeli dejanskega sistema. V modelu so bili tokovi in napetosti manjše od tistega v resnici. Danes se analogne modele mreže še vedno uporablja, saj po svoji strukturi zagotavljajo delovanje v realnem času, poleg tega pa tudi pomanjšani tokovi in napetosti lahko zagotavljajo za preizkušano opremo v zaprti zanki dovolj realne delovne pogoje, saj se v pravem elektroenergetskem sistemu uporablja tokovne in napetostne zaščitne transformatorje. Njihova slabost je v tem, da je potrebno imeti fizični model in veliko znanja o stvareh, ki s problemom, ki ga opazujemo niso neposredno povezane. Tudi strošek TNA je lahko zelo visok, saj morajo biti elementi kalibrirani. Pozneje so pogosto uporabljali hibridne simulacije, kjer so fizični modeli mreže, kot je prej omenjen analogni TNA nadomeščali en del omrežja, dopolnitev k temu modelu je bila pa simulirana. Ta način se imenuje hibridna simulacija in je bolj ali manj opuščen. Pozneje se je z razvojem mikroprocesorjev in predvsem signalnih procesorjev s tehnologijo plavajoče vejice (floating point) začelo pojavljati digitalno dinamično simuliranje. Prve digitalne simulacije v realnem času so se izvajale v osemdesetih letih prejšnjega stoletja na specializiranih digitalnih simulatorjih, ki so bili narejeni posebej za ta namen. Uporabljali so jih za preizkuse opreme v zaprti zanki. Pozneje so se začela pojavljati programska orodja za superračunalnike. Z razvojem hitrejših osebnih računalnikov , predvsem pa z večjedrnimi procesorji so postale simulacije v realnem času veliko bolj izvedljive na serijsko izdelanih osebnih računalnikih in simulatorjih. Zasnova simulacije na osebnih računalnikih se v literaturi [1] običajno omenja kot COTS (Commercial Off The Shelf). Danes obstaja več programskih orodij za simuliranje v Pomurska obzorja 1/ 2014/ 2 | 67 Rafael MIHALIČ: SIMULIRANJE ELEKTROMAGNETSKEGA SISTEMA – ZAKAJ IN KAKO realnem času in posledično opravljanje zaprtozančnih preizkusov, ki izvajajo izračune simulacije na osebnem računalniku. Za njihovo izvedbo je običajno potrebno imeti nek poseben vmesnik, ki ojači signal, zagotavlja A/D pretvorbo in prevzame nase breme komuniciranja in prevajanja podatkov. Odločilen pospešek razvoju tehnike, ki omogoča simulacijo v realnem času je dala potreba po treniranju vojaških pilotov in pilotov potniških letal. V taki simulaciji tečejo podatki do osebe, ki se ukvarja s simulatorjem in nazaj, zato je tudi to simulacija v zaprti zanki. Sicer se naloge pilota in sistemskega operaterja elektroenergetskega sistema skoraj v celoti razlikujejo, vendar so odločitve v izrednih razmerah obeh zasnovane na preteklih izkušnjah, na priučenih reakcijah iz treninga in na razumevanju dogajanja. V preteklosti so že izvajali projekte, kjer so uporabljali simulacije v realnem času tudi za treniranje sistemskih operaterjev [2]. Kot je bilo že omenjeno je ključni pogoj za izvedbo simulacije v realnem času, da je izvajanje matematičnih operacij dovolj hitro. Pri tem imamo načeloma 2 možnosti, in sicer, da bodisi izvedemo računske operacije zelo hitro ali pa uporabimo paralelno procesiranje. Prvi pristop ima slabost, da potrebujemo zelo hitre računalnike (superračunalnike), in da sistem ne more biti poljubno velik. V drugem primeru pa seveda potrebujemo še vedno hitre računalnike, vendar lahko sistem širimo z dodajanjem novih procesorskih in komunikacijskih enot. Vendar pa je treba uporabiti metode, ki paralelno procesiranje omogočajo. 2.1. Matematična podlaga paralelnega procesiranje pri simulaciji EES Računalniško zmogljivost paralelnega procesiranja je moč izkoristiti le, če reševanje problema na nek način porazdelimo med več procesorjev. Ključen element je t.i. homogen vod s porazdeljenimi parametri, zato si glavne ideje oglejmo nekoliko pobliže. Shemo infinitezimalno majhnega delčka voda in pripadajočih električnih veličin prikazuje slika 2. Za električni vod s porazdeljenimi parametri brez izgub veljata diferencialni enačbi: Ob vpeljavi : n= 1 LC u = i = 1 lc = 1 (hitrost širjenja) in (3) (valovna impedanca) (4) k m0 e0 l = Z0 c je mogoče s parcialnim odvajanjem posamezne enačbe (1) oz. (2) in s preureditvijo dobiti d'Alembertovo enačbo drugega reda. Nadalje s substitucijo spremenljivk, odvajanjem in vstavitvijo odvodov v prvotno d'Alembertovo enačbo dobimo pogoj, da mora biti mešani odvod toka oz. napetosti po novih spremenljivkah enak nič. To pomeni, da je analitična rešitev za tok oz. napetost seštevek nekih ločenih funkcij vpeljanih spremenljivk. Funkciji sta lahko poljubni, pomembno je le, da sta za enake vrednosti spremenljivk enaki. Dobimo: i ( x, t ) = f1 ( x - vt ) + f2 ( x + vt ) . (5) u ( x, t ) = Z 0 f1 ( x - vt ) - Z 0 f2 ( x + vt ) . (6) Pri tem velja: u napetost i tok t čas x krajevna spremenljivka (lokacija), k faktor ( permeabilnosti in dielektričnosti npr. kabel), μ0 permeabilnost praznega prostora, ε0 dielektričnost praznega prostora, v hitrost širjenja valov, l induktivnost na 1 km voda (L/ dolžina v km), c kapacitivnost na 1 km voda (C/ dolžina v km), f1,f2 funkciji toka, Z0 karakteristična impedanca voda. u u1(x1,t1)=u2(x2,t2) t1 t2 x1 x2 (1) (2) x Slika 3. Potujoč napetostni val. i(x,t) i(x+dx,t) l dx + + u(x,t) c dx - u(x+dx,t) - Slika 2. Diferencialni odsek brezizgubega voda 68 | Pomurska obzorja 1/ 2014/ 2 V bistvu dobimo potujoč val napetosti ali toka, pri čemer en člen predstavlja val v eno, drug člen pa val v drugo smer (v enačbah 5, 6 – t.j. "desno oz. levo"). Za posamično vozlišče (recimo A oz B) z upoštevanjem potovalnega časa voda lahko zapišemo enačbi toka in napetosti v vozliščih in dobimo: iA (t ) = 1 uA (t ) + iB (t - t ) Z0 (7) iB (t ) = 1 uB (t ) + iA (t - t ) . Z0 (8) Rafael MIHALIČ: SIMULIRANJE ELEKTROMAGNETSKEGA SISTEMA – ZAKAJ IN KAKO  čas potovanja vala Zadnji dve enačbi že vsebujeta mehanizem za razdelitev celotnega sistema na paralelne podsisteme. Vsaka od obeh enačb povezuje sedanji tok v enem vozlišču s sedanjo napetostjo v tem vozlišču in preteklo napetostjo oz. preteklim tokom v nasprotnem vozlišču voda s porazdeljenimi parametri. Zato se lahko računa pri takem vodu sistem ločeno na vsaki strani voda. Rešitev toka in napetosti na eni strani voda v nekem trenutku namreč ne vpliva na rešitev toka in napetosti na drugi strani voda v tem istem trenutku. Odvisnost med vozlišči je zajeta v delih, ki so odvisni od časa zamaknjenega za potovalni čas v preteklost. Vsak procesor lahko torej obravnava del sistema neodvisno, vendar pa mora pred naslednjim korakom izračuna dobiti podatek o razmerah na drugi strani voda v pravkar minulem časovnem koraku in nato izračunati razmere v naslednjem koraku. V realnem času se mora torej izvršiti, ne samo izračun, pač pa tudi komunikacija. Izračune izvajajo hitri RISC procesorji, komunikacija pa poteka preko optičnih vlaken. Pri do sedaj doseženih hitrostih izračunov in komunikacije (računski korak 50 s) se lahko ločijo deli omrežij, ki jih povezujejo vodi 15 km (in več). 3. Računalnik za simulacijo dinamike elektroenergetskih sistemov v realnem času Na svetu je več izdelovalcev strojne opreme za simulacijo EES v realnem času. RDTS Technologies iz Kanade [5] je med njimi vodilna v proizvodnji kakovostnih simulatorjev, ki so modularne izvedbe in izjemnih karakteristik ter referenc povsod po svetu. Simulatorji RTDS so sestavljeni iz enot, poimenovanih kot »Rack« (polica). Vsaka polica vsebuje lahko do šest procesorskih kartic (Giga Processor Card– GPC ali najnovejša serija kartic PB5) in vsaj eno kartico GTWIF (angl. Giga Transceiver Workstation Inter Facecard). Slednja je nujno potrebna za komunikacijo med osebnim računalnikom in simulatorjem, ki jo izvaja prek komunikacije TCP/IP (Ethernet). Simulacija v realnem času je lahko kakovostna le s čim manjšim časovnim korakom integracije, pri čemer simulator RTDS lahko doseže časovni korak integracije tudi 50 μs (za tako imenovani »Small-time-step modeling« tudi do 1.5 μs). S takšno hitrostjo izračuna pa je nemogoče pričakovati, da bo osebni računalnik prek povezave Ethernet sproti prebiral vse podatke iz simulatorja in jih tudi v celoti v realnem času prikazoval na zaslon. Zato se prek modula RSCAD/RunTime opazujejo spremenljivke v celoti zgolj na ukaz. Simulator, ki ga imamo na voljo na Fakulteti za elektrotehniko v Ljubljani, vsebuje eno polico, na kateri so štiri procesorske kartice (vsaka s po dvema procesorjema) in ena kartica GTWIF. Poleg strojne opreme je sestavni del simulatorja tudi programska oprema RSCAD, prek katere je mogoče samo strojno opremo upravljati, ter je nameščena na osebnem računalniku. Simulator in osebni računalnik sta lahko povezana neposredno ali prek računalniškega omrežja, če simulator uporablja več uporabnikov. Spremljanje časovnega spreminjanja veličin v modelu EES, ki ga v realnem času izvaja simulator RTDS, je mogoč na več načinov. Prvi je prek modula RSCAD/RunTime (na primer prek grafa), drugi pa prek opazovanja električnega signala (npr. z osciloskopom), kar je pravzaprav posebnost takšnega simulatorja. V ta namen je mogoča uporaba analogne (angl. GT Analogue Outputcard - GTAO) in digitalne izhodne kartice (angl. GT Digital Outputcard - GTDO) ali pa vmesnika na prednji plošči procesorske kartice GPC (slika 4a) [4]. Slednji je namenjen predvsem opazovanju signalov z osciloskopom, medtem ko sta izhodni kartici sposobni proizvesti visokoresolucijski signal, potreben predvsem za povezovanje simulatorja in krmilnih oziroma zaščitnih naprav. a) b) Slika 4. a) kartica GPC z vmesnikom na prednji plošči, b) ojačevalec (ojačevalnik signala) Kot primer smo izvedli vklop dodatnega voda vzporedno k že obratujočemu vodu, pri čemer pa smo sočasno merili tok faze A skozi vklopljeni vod tako v modulu RSCAD/RunTime (slika 5- zgoraj) kot tudi z osciloskopom prek vmesnika na prednji plošči procesorske kartice GPC (slika 5- spodaj). Kot že zapisano, je slednja v primerjavi skartico GTAO sposobna proizvesti signal s slabšim dinamičnim sledenjem dejanskemu signalu [5], kar je razvidno tudi iz primerjave obeh grafov – slika 5. Treba se je zavedati, da je razlika med grafoma nastala kljub uporabljenemu ročnemu osciloskopu s frekvenco vzorčenja 40 MHz. Iz tega sledi, da je treba, če je potreben izvoz signala z boljšim dinamičnim sledenjem, uporabiti kartico GTAO. I (kA) Pri tem velja: Čas (sek) 85 ms Slika 5. Časovni poteki toka faze A po vklopljenem vodu vzporedno k že obratujočemu vodu –RSCAD/RunTime (zgoraj) in ročni osciloskop (spodaj). Pomurska obzorja 1/ 2014/ 2 | 69 Rafael MIHALIČ: SIMULIRANJE ELEKTROMAGNETSKEGA SISTEMA – ZAKAJ IN KAKO 3.1. Simulacija v zaprti zanki z zaščitnim relejem Kot rečeno je posebnost RTDS simulatorja možnost analize/testiranja določenih naprav v zaprti zanki z modelom preostalega EES (angl. Closed Loop Testing). Pri takšnih simulacijah se del modela EES nadomesti z dejansko napravo, ki se priključi na simulator RTDS. To mu med drugim omogočata analogna in digitalna vhodna kartica. Poleg teh je mogoče binarne signale v simulator uvesti tudi prek vmesniškega panela na prednji plošči simulatorja, ki je namenjen predvsem prenosu izklopnih signalov zaščitne opreme. Kartica GTAO ima razpon električnega signala na svojih sponkah v mejah ±10 V, pri čemer pa je njena tokovna sposobnost do 25 mA. Iz tega sledi, da je pri testiranju zaščitnega nadtokovnega releja potrebna uporaba tokovnega ojačevalnika. Za demonstracijo smo uporabili EES je sestavljen iz dveh sinhronskih generatorjev, ki prek štirih daljnovodov napajata dve odjemni zbiralki. Uporabljeni nadtokovni rele, priključen na sekundarne sponke modela tokovnega transformatorja (CT) varuje vod L2 ter v primeru delovanja odklopi oba odklopnika BR_1 in BR_2. 1 L1 Load_A REF630 L2 BR_1 Gen_1 Kot že rečeno ima kartica GTAO razpon električnega signala na svojih sponkah v mejah ±10 V. Ta signal je potrebno prilagoditi nivoju vhodnega signala naprave. V konkretnem primeru gre za nadtokovni rele, torej je vhodni signal tokovni, reda nekaj amperov. Tukaj trčimo na problem prilagoditve signala z ojačevalnikom. Moč je sicer uporabiti standardni generator signalov za testiranje zaščitnih relejev (npr. Omicron CMC 356), ki ga ustrezno krmilimo [7], vendar se moramo ob tem zavedati, da pride do zakasnitve signala reda 100 ms. Za relativno počasne pojave je to še sprejemljivo, za analizo hitrih pojavov in testiranje krmilnikov močnostne elektronike pa nikakor ne. Za verno ponazoritev hitrih pojavov nujno potrebujemo ustrezen ojačevalnik, katerega zakasnitev se nahaja v področju reda s, torej je reda 1000 x hitrejši. Ob tem ne smemo pozabiti, da vsak 3-fazni signal zahteva 3 vhode in 3 izhode. V primeru npr. distančnega releja to pomeni 3 tokovne in 3 napetostne izhode (in ravno toliko vhodov) ojačevalnika. Rezultat simulacije kratkega stika in krmiljenja izklopa z dejansko napravo prikazuje slika 8. BR_2 Gen_2 KS 2 3 L3 L4 4 Load_B Slika 6. Model testnega EES. V nekem trenutku, ki ga interno sprožimo prek modula RSCAD/RunTime, pride na zbiralki 1 do tripolnega kratkega stika. Rele, ki ima najenostavnejšo tokovno neodvisno izklopno karakteristiko (slika 7a) prejme signal kratkostičnega toka in po nastavljenem času odda izklopni signal. Ta signal se uvede nazaj v računalnik in stikali BR_1 in BR_2 v modelu se ob prehodu toka skozi nič (simulacija obloka) izklopita. Shema pretoka signalov je ponazorjena na sliki 7b. 4. Sklep Cilj predstavljenega dela je predvsem ilustrirati pomen simulacije EES in demonstrirati uporabnost nove velike pridobitve laboratorijev LPEE in LEON na Katedri za elektroenergetske sisteme in naprave, Fakultete za elektrotehniko, UL, t.j. simulatorja za digitalno simulacijo elektroenergetskega sistema v realnem času RTDS. Slednje smo b) REF 630 a) Slika 8. Časovni potek tokov na sekundarju tokovnika. TRIP IA IB IC GND VA VB VC GND CMC 356 GTAO IA IB IC CT Model EES OPEN/ CLOSE FPI RTDS Simulator Slika 7. a) Nastavitev nadtokovne karakteristike zaščitnega releja REF630 izdelovalca ABB b) Pretok signalov testnega sistema 70 | Pomurska obzorja 1/ 2014/ 2 Rafael MIHALIČ: SIMULIRANJE ELEKTROMAGNETSKEGA SISTEMA – ZAKAJ IN KAKO izvedli na primeru preizkušanja opreme, preden bi jo vključili v dejanski elektroenergetski sistem. Tak simulator je v svetu nepogrešljivo orodje pri analizi obnašanja krmilne opreme (na primer regulacija pretvornikov HVDC) preden jo vgradijo v EES. V realnem sistemu namreč tovrstnih preizkusov ni mogoče izvesti. Kot razlog za to trditev si predstavljajmo na primer regulator 5000 MW 800 kV aplikacije HVDC, ki ga je treba preizkusiti za primere okvar v EES zelo visokih napetosti [8]. Simulator RTDS omogoča testiranje dejanske opreme v zaprti zanki z digitalnim modelom preostalega elektroenergetskega sistema. Časovni korak izračuna elektromagnetnih pojavov, ki jih simulator RTDS omogoča, se giblje nekje do 50 μs, medtem ko je mogoče za modeliranje elementov močnostne elektronike (na primer tiristorski pretvorniki in naprave FACTS) ta korak znižati tudi do1.5 μs. Ker lahko vsak procesor obdeluje le določeno število elementov, se uporabi Dommelova metoda modeliranja homogenih vodov [9] in na ta način sistem porazdeli na območja, ki jih lahko v okviru enega integracijskega koraka obravnavamo neodvisno od ostalega sistema, če imamo na razpolago podatke na zunanjih mejah omenjenega območja iz preteklosti (zakasnitev za potovalni čas svetlobe vzdolž povezav tega območja z ostalim sistemom). Kot praktični primer uporabe smo izvedli preizkus enostavne nadtokovne funkcije releja REF 630 izdelovalca ABB. Prednost testiranja relejev prek simulatorja RTDS je v tem, da se preizkusi rele v razmerah, ki so krepko bliže dejanskim razmeram v EES, skupaj s harmonskim popačenjem vhodnih signalov, nihanji frekvence in podobnim. Pri preizkusu se je izkazalo, da je v večini primerov poleg samega simulatorja potrebna tudi uporaba hitrega ojačevalnika tokovnih in napetostnih signalov, kar lahko predstavlja nezanemarljiv problem. 5. Ankush Saran, Padmavathy Kankanala, Anurag K. Srivastava, Noel N. Schulz, Designing and Testing Protective Overcurrent Relay Using Real Time Digital Simulation, RTDS Technologies internal library. 6. OMICRON, CMC356 Technicaldata, dostopno na http://www.omicron.at/fileadmin/user_upload/files/pdf/en/ CMC-356-Brochure-ENU.pdfna dan 09.02.2012. 7. Omicron: Programming interface for CMC test systems, CMENG.AE.6, 2007 8. H. Duchen, M. Lagerkvist, R. Kuffel, R.P. Wierckx, HVDC Simulation and Control System Testing Using a Real-Time Digital Simulator (RTDS), RTDS Technologies internal library. 9. H. W. Dommel, Digital Computer Solution of Electromagnetic Transients in Single and Multiphase Networks, IEEE Transactions on Power Apparatus and Systems, vol. PAS-88, N0. 4, April 1969, pp. 388-399. Zahvala To delo je financirala Javna agencija za raziskovalno dejavnost Republike Slovenije (ARRS) v sklopu programske skupine Elektroenergetski sistemi, P2-356. Hvala! Literatura 1. J. Bélanger, P. Venne, J.N. Paquin , The What, Where and Why of Real-Time Simulation, publikacija podjetja OpalRT dostopna na spletnem naslovu : http://www.opalrt.com/sites/default/files/ technical_papers/PES-GMTutorial_04%20-%20Real%20Time%20Simulation.pdf 2. U. Spaniel, M. Kreutz, C Roggatz ; Simulator based operator training-ensuring the quality of power system operation; IEEE Power System Management and Control, 2002. Fifth International Conference on (Conf. Publ. No. 488) 3. RTDS Technologies, Real time digital simulation for the power industry - manual set, RSCAD version 2.024.2. 4. RTDS Technologies, /index.html. http://www.rtds.com/index Pomurska obzorja 1/ 2014/ 2 | 71 Tehnika Zoran Žunič* Mešana metoda robnih in končnih elementov za reševanje hitrostno vrtinčne formulacije Navier-Stokesovih enačb POVZETEK V prispevku je predstavljena numerična metoda za reševanje nestisljivega toka viskozne newtonske tekočine z uporabo hitrostnovrtinčne formulacije ohranitvenih enačb. Izpeljali smo robno integralsko enačbo kinematike v tangentni obliki, s čemer je bil dosežen natančnejši izračun robnih vrednosti vrtinčnosti. Po metodi končnih elementov smo izpeljali še integralsko enačbo kinematike, s katero smo izračunali neznane območne hitrosti in pa enačbo kinetike, s katero smo izračunali neznane območne vrtinčnosti. Vse dobljene sisteme diskretnih enačb smo povezali v računski algoritem, ki z metodo Picardove iteracije nesklopljeno reši vezan sistem nelinearnih enačb. Da bi pokazali robustnost in uporabnost predstavljene metode in istočasno izločili možnost napak v programu, smo uporabili primer gnanega toka v prostorski kotanji. Ključne besede: prostorski tok viskozne tekočine, Navier-Stokesove enačbe, hitrostno-vrtinčna formulacija, metoda robnih elementov, metoda končnih elementov. 1. Uvod Pri razvoju na področju numeričnih metod za obravnavanje tokov viskoznih tekočin se običajno osredotočamo na aproksimacijske metode za reševanje Navier-Stokesovih enačb. Večina teh pristopov uporablja le eno od aproksimacijskih metod, na primer metodo končnik razlik, metodo končnih prostornin, metodo končnih elementov ali metodo robnih elementov. Vsaka od teh metod ima svoje prednosti in slabosti, zato se je začel razvoj tudi na področju istočasne uporabe večih metod za reševanje istega sistema enačb. Z uporabo hitrostno-vrtinčne formulacija so se v preteklosti ukvarjali že Liu [1], ki je uporabil metodo končnih razlik, Guevremont in dr. [2], [3] ter Wong in Baker [4], ki so uporabili metodo končnih elementov. Na področju robnih elementov so delovali Žagar in Škerget [5], Hriberšek in Škerget [6], Ramšak in dr. [7] ter Ravnik in dr. [8]. Nekaj raziskovalcev se lotilo tudi uporabe večih metod, med katere spadajo Young in dr. [9] ki so uporabljali metodo robnih elementov in metodo končnih elementov vendar s formulacijo z osnovnimi spremenljivkami, ter Brown in Ingber [10], ki sta notranje hitrosti računala z posplošeno Helmhotzovo metodo dekompozicije, transport vrtinčnosti pa z Galerkinovo metodo končnih elementov. Prikazano delo je kombinacija reševanja vektorske Poissonove hitrostne enačbe z metodo robnih elementov za izračun robnih vrtinčnosti in uporabe metode končnih elementov za reševanje vektorske Poissonove hitrostne enačbe za notranje hitrosti in transportne enačbe vrtinčnosti. Prvi pristop omogoča * AVL – AST d.o.o. / Ulica kneza Koclja 22, SI-2000 Maribor E-naslov: zoran.zunic@amis.net Dr. Zoran Žunič, Lorbekova 15, 2341 Limbuš, Slovenia; Tel.: +386-2-250-85-16 72 | Pomurska obzorja 1/ 2014/ 2 implicitni izračun robnih vrtinčnosti, medtem ko drugi pristop reši težavo s prevelikimi zahtevami po računalniškem pomnilniku. 2. Hitrostno vrtinčna formulacija Tok gibanja zvezne tekočine opisujeta enačbi ohranitve gibalne količine in ohranitve mase, ki sta zapisani z osnovnimi spremenljivkami (tlaki in hitrostmi). Pri reševanju problemov nestisljivih tokov, pri katerih predpostavimo, da je gostota konstantna, se pojavi težava zaradi odsotnost tlaka in gostote v enačbi ohranitve mase. Težavo rešimo z vpeljavo pojma vrtinčnosti,      (1)     v,    0 s pomočjo katerega preoblikujemo osnovne enačbe v hitrostno – vrtinčno obliko (Škerget in dr. [7], Ravnik in dr. [8]). Dobimo Poissonovo enačbo za hitrost (enačbo kinematike)     2v      0 , (2) in ob predpostavki Newtonske tekočine še transportno enačbo za vrtinčnost (enačbo kinetike)           (v  )  (  )v   2 . t (3) Iščemo rešitev enačb (2) in (3), ki zadovoljuje začetne in robne pogoje v računskem območju. 3. Metoda robnih elementov Enačbo (2) zapišemo kot singularno robno integralsko enačbo za vektor hitrosti s pomočjo Greenovih teoremov (Wrobel [11])          c ( )v ( )   v ( n   )u * d    u * ( n   )v d    (   )u * d (4)    Zoran ŽUNIČ: MEŠANA METODA ROBNIH IN KONČNIH ELEMENTOV ZA REŠEVANJE… kjer u predstavlja eliptično Laplaceovo osnovno rešitev,  je 3. Zanka zaradi nelinearnosti. izvorna točka na robu računskega območja, c je geometrijski  koeficient in n normala na rob računskega območja. Enačbo (4) z Gausovim divergenčnim teoremom pretvorimo v tangentno obliko, s čimer se znebimo odvodov hitrostnega in vrtinčnega polja (Škerget in dr. [7], Ravnik in dr. [8]), ki jih prenesemo na osnovno rešitev          c ( )v ( )   (u *  n )v d    (u *  n )  v d      u * d 4. Kinematika toka: *   a. izračunamo robne vrtinčnosti (5) z metodo robnih elementov, b. izračunamo notranje hitrosti (6) z metodo končnih elementov. 5. Kinetika toka, transport vrtinčnosti: a. izračunamo notranje vrtinčnosti (7) z metodo končnih elementov, b. podrelaksiramo rešitev.  (5) 6. Preverimo napako rešitve, če je večja od tolerance gremo na korak 3. 4. Metoda končnih elementov Z uporabo Galerkinove metode utežnih ostankov in uporabo Greenovega prvega teorema v enačbi (2) dobimo šibko obliko Poissonove enačbe za hitrostni vector           (N   )v d   N ( n   )v d    N (   ) d  0 (6)    kjer nam N predstavlja utežno funkcijo. Enak postopek ponovimo z enačbo (3), da dobimo šibko obliko prenosne enačbe za vrtinčnost        1 N d   N (v   )  d   N (   ) v d   t           1   (N   ) d    N (n   ) d    N 1 d  t    (7) kjer je  t velikost časovnega koraka in podpis   1 vrednost v prejšnjem časovnem koraku. Velikosti integralov potrebnih za reševanje enačb (6) in (7) so veliko manjše od velikosti integralov v enačbi (5), kar nam omogoča hitrejše reševanje problemov pri enaki gostoti računske mreže, oziroma obravnavanje večjih (kompleksnejših) problemov. 5. Računski postopek Računski postopek je potekal v naslednjem vrstnem redu: 1. Izberemo začetno hitrostno polje in izračunamo začetno vrtinčno polje z enačbo (1). Začetni čas =0. 7. Konec časovnega koraka: a. shranimo vrednosti hitrosti in vrtinčnosti, b. če je čas manjši od predpisanega gremo na korak 2. 8. Konec izračuna. 6. Primerjava rezultatov Gnan tok v prostorski kotanji je eden najpogosteje uporabljanih testnih primerov za testiranje laminarnih tokov viskozne Newtonske tekočine. Kljub enostavni geometriji, se v njem pojavljajo vsi fenomeni kompeksnih tokov, kot so vrtinci, sekundarni vrtinci, prostorski vzorci, kaotično gibanje, nestabilnosti in turbulence (Shankar in Despande [12]). V našem primeru smo uporabili kotanjo v obliki kocke s stranicami dolžine 1 enote. Reynoldsovo število, ki je v tem primeru edino karakteristično število, je bilo zasnovano na dolžini stranice kocke in hitrosti vlečenja zgornje ploskve (Re = vx L/) in je znašalo Re=100 in 1000. Velikost uporabljenih mrež je znašala 8 3, 163 in 323 elementov. Vse mreže so bile zgoščene proti zunanjim stenam s faktorjem 8 (slika 1a). Predpisali smo naslednje robne pogoje:  zgornja ploskev (𝑧 = 1): 𝑣𝑥 = 1, 𝑣𝑦 = 𝑣𝑧 = 0,  ostale ploskve: 𝑣𝑥 = 𝑣𝑦 = 𝑣𝑧 = 0. Vse testne primere smo izračunali z začetnimi pogoji 𝑣𝑥 = 𝑣𝑦 = 𝑣𝑧 = 0. 2. Časovna zanka. Slika 1. a) Računska mreža, b) Hitrostni profili za Re=1000. Pomurska obzorja 1/ 2014/ 2 | 73 Zoran ŽUNIČ: MEŠANA METODA ROBNIH IN KONČNIH ELEMENTOV ZA REŠEVANJE… Na sliki 1b je prikazana primerjava profilov hitrosti, dobljenih na različno gostih mrežah, z rezultati Yanga in dr. [13]. Primerjali smo hitrostne profile vx in vy po navpični in vodoravni srednjici kocke. Analiza z zgoščevanjem mreže je pokazala, da je najredkejša mreža pregroba, medtem ko se gostejši mreži dobro ujemata z referenčnimi.  Slika 2 prikazuje izopovršine absolutne hitrosti | v |=0.13 za različne vrednosti Reynoldsovega števila. Vidimo, da se hitro vrteče se jedro toka pomika proti zunanjim stenam, obenem pa se debelina hitrega pasu zmanjšuje. Slika 3 prikazuje prikazuje vrtinčno polje na različnih presekih kotanje (x = 0.5, y = 0.5 in z = 0.5) za vrednost Reynoldsovega števila Re = 1000. Iz izolinij vrtinčnosti je razvidno, da je tokovno polje pri izbrani vrednosti Reynoldsovega števila simetrično preko xz ravnine.  Slika 2. Izopovršine hitrosti | v |=0.13. Levo Re=100, desno Re=1000. Slika 3. Izolinije vrtinčnosti za Re=1000. Zgoraj levo x, zgoraj desno y, spodaj z . 74 | Pomurska obzorja 1/ 2014/ 2 Zoran ŽUNIČ: MEŠANA METODA ROBNIH IN KONČNIH ELEMENTOV ZA REŠEVANJE… 7. Zaključki V prispevku smo prikazali postopek reševanja NavierStokesovih enačb zapisanih s hitrostno-vrtinčno formulacijo. Robne vrtinčnosti smo izračunali implicitno z uporabo metode robnih elementov, medtem ko smo notranje hitrosti in notranje vrtinčnosti izračunali z uporabo metode končnih elementov. Na tak način smo ohranili visoko natačnost izračuna robnih vrtinčnosti, prihranili pa smo pri količini potrebnega računalniškega pomnilnika. Za preverjanje pravilnosti in prikaz natančnosti uporabljene metode smo uporabili primer gnanega toka v prostorski kotanji. Primerjava dobljenih rezultatov z rezultati iz literature je pokazala odlično ujemanje. 12. P. N. Shankar, M. D. Deshpande, Fluid mechanics in the driven cavity, Annu. Rev. Fluid Mech. 32, 93–136, 2000. 13. J. Y. Yang, S. C. Yang, Y. N. Chen, C. A. Hsu, Implicit weighted ENO schemes for the three-dimensional incompressible Navier-Stokes equations, Journal of Computational Physics 146, 464–487, 1998. Literatura 1. C. H. Liu, Numerical solution of three-dimensional NavierStokes equations by a velocity-vorticity method, Int. J. Numer. Meth. Fluids 35, 533–557, 2001. 2. G. Guevremont, W. G. Habashi, Finite element solution of the Navier-Stokes equations by a velocity-vorticity method, Int. J. Numer. Meth. Fluids 10, 461–475, 1990. 3. G. Guevremont, W. G. Habashi, P. L. Kotiuga, M. M. Hafez, Finite element solution of the 3D compressible Navier-Stokes equations by a velocity-vorticity method, Journal of Computational Physics 107 (1), 176–187, 1993. 4. K. L. Wong, A. J. Baker, A 3D incompressible NavierStokes velocity-vorticity weak form finite element algorithm, Int. J. Numer. Meth. Fluids 38, 99–123, 2002. 5. I. Žagar, L. Škerget, Boundary elemens for time dependent 3-D laminar viscous fluid flow, Journal of mechanical engineering 35 (10/12), 160–163, 1989. 6. M. Hriberšek, L. Škerget, Iterative methods in solving Navier-Stokes equations by the boundary element method, Int. J. Numer. Methods Engnrg 39, 115–139, 1996. 7. M. Ramšak, L. Škerget, M. Hriberšek, Z. Žunič, A multidomain boundary element method for unsteady laminar flow using stream function – vorticity equations, Eng. Anal. Bound. Elem. 29, 1–14, 2005. 8. J. Ravnik, L. Škerget, M. Hriberšek, Two-dimensional velocity-vorticity based LES for the solution of natural convection in a differentially heated enclosure by wavelet transform based BEM and FEM, Eng. Anal. Bound. Elem. 30, 671–686, 2006. 9. D. L. Young, J. L. Huang, T. I. Eldho, Simulation of laminar vortex shedding flow past cylinders using a coupled BEM and FEM model, Comput. Methods Appl. Mech. Engrg. 190, 5975–5998, 2001. 10. M. J. Brown, M. S. Ingber, Parallelization of a vorticity formulation for the analysis of incompressible viscous fluid flows, Int. J. Numer. Meth. Fluids 39, 979–999, 2002. 11. L. C. Wrobel, The Boundary Element Method, John Willey & Sons, LTD, 2002. Pomurska obzorja 1/ 2014/ 2 | 75 Tehnika Iztok Fister1,*, Janez Brest1, Iztok Fister Jr.1, Karin Ljubič2 Algoritmi po vzorih iz narave POVZETEK Narava je bila vedno vzor znanstvenikom pri reševanju najtežjih problemov v računalništvu, matematiki, medicini in industriji. Vzori iz narave ponujajo nove smernice za razvoj računalniških algoritmov. Ti algoritmi lahko temeljijo na obnašanju inteligence rojev, oz. drugih bioloških sistemov ali na simulaciji fizikalnih ali kemijskih procesov. V tem članku predstavljamo značilnosti vsakega od treh omenjenih vzorov, ki so vplivali na nastanek novih računalniških algoritmov, omenimo njihove avtorje in članke, v katerih so predstavili svoja dela. S tem pregledom želimo pomagati razvijalcem programske opreme pri lažji izbiri orodij za reševanje njihovih problemov iz prakse. Ključne besede: algoritem, inteligenca rojev, biološki sistemi, kresnice, netopirji. 1. Uvod Reševanje optimizacijskih problemov je težka in zelo zapletena naloga, saj je večina teh problemov NP-težkih (Garey & Johnson, 1979), kar pomeni, da se časovna zahtevnost reševanja povečuje eksponentno z večanjem naloge. To za uporabnika pomeni, da lahko pričakuje rezultate optimizacije šele po dolgotrajnem računanju, v najslabšem primeru pa to dejansko pomeni, da rešitve nikoli ne najdemo. Ker je izčrpno iskanje (pregled vseh možnih rešitev v prostoru preiskovanja) časovno preveč zahtevno, poskušamo najti optimalne rešitve v realnem času s pomočjo aproksimativnih metod »generiraj in preveri« (angl. generate and test). V zgodovini je bilo razvitih že veliko uspešnih mehanizmov za reševanje teh problemov, od različnih matematičnih modelov, linearnega programiranja, ipd… V sodobnem času, ko postajajo računalniki vse zmogljivejši, raziskovalci pri razvoju novih vrst algoritmov iščejo rešitve iz vzorov iz narave. Pri tem poskušajo na računalniku posnemati različne procese v naravi in te uporabiti za reševanje težkih problemov v realnem svetu. Najstarejši predstavniki te vrste so evolucijski algoritmi (Eiben & Smith, 2003). Ti algoritmi posnemajo biološko teorijo boja za obstanek , ki jo je razvil Charles Darwin (1859). V devetdesetih letih prejšnjega stoletja se je tem pridružila družina algoritmov, ki temeljijo na kolektivnem obnašanju v samoorganiziranih in decentraliziranih družinah žuželk oz. živali, kot npr. mravelj, termitov, čebel, jat ptičev, ipd… Vzore za reševanje težkih problemov iščemo tudi v drugih naravoslovnih znanostih, kot npr. kemiji in fiziki, v nekaterih družboslovnih znanostih in celo v glasbi. 1 Univerza v Mariboru, Fakulteta za elektrotehniko, računalništvo in informatiko, Smetanova 17, 2000 Maribor 2Univerza v Mariboru, Medicinska fakulteta, Taborska 8, 2000 Maribor E-naslovi: iztok.fister@uni-mb.si; janez.brest@uni-mb.si; iztok.fister2@uni-mb.si; karin.ljubicc@gmail.com 76 | Pomurska obzorja 1/ 2014/ 2 Nekateri od teh algoritmov so se v tem času že uveljavili, saj smo jih uspešno uporabili pri reševanju težkih problemov v realnem svetu (npr. optimizacija z delci rojev, algoritmi na osnovi obnašanja kresnic, netopirjev, kukavic, ipd…), drugi se šele uveljavljajo, nekateri pa ostajajo bolj ali manj na nivoju ideje. Za uveljavitev novega algoritma niso pomembni samo dobri rezultati, ki jih ta dosega pri reševanju določenega problema, ampak je potrebno zanje, prek konferenc oz. objav v znanstvenih revijah, zbuditi pozornost mednarodne strokovne javnosti. V tem članku smo zbrali prek 40 vrst algoritmov iz različnih področij uporabe, ki delujejo po vzorih iz narave. Te algoritme poskušamo klasificirati v različne kategorije glede na njihov izvor in ugotoviti, na katerih področjih so se še posebej uveljavili. Članek temelji na izhodiščih, ki so jih postavili Fister Jr. in drugi (2013). Struktura članka v nadaljevanju je naslednja: v poglavju 2 predstavimo vzore, ki navdihujejo znanstvenike pri razvoju novih vrst algoritmov. Pri tem se še posebej osredotočamo na algoritme po vzorih iz narave. Poglavje 3 se ukvarja s klasifikacijo algoritmov po vzorih iz narave. Članek zaključimo s povzetkom opravljenega dela in začrtamo smernice za nadaljnji razvoj. 2. Vzori za razvoj novih vrst algoritmov Narava je vedno navdihovala računalniške strokovnjake pri razvoju novih algoritmov, zato ni čudno, da največ vzorov za razvoj novih vrst algoritmov prihaja ravno iz biologije. Najpomembnejši algoritmi te vrste so evolucijski algoritmi (Eiben in Smith 2003), katerih vzor je Darwinov boj za obstanek (1859), kjer imajo v tekmi za preživetje več možnosti najuspešnejši posamezniki, ki se znajo prilagoditi spremembam v okolju. Iztok FISTER et al.: ALGORITMI PO VZORIH IZ NARAVE Tudi inteligenca rojev je nastala po vzorih iz biologije (Blum in Merkle 2008). V tem primeru predstavljajo njihov vzor na videz enostavna bitja sposobna izvajanja prirojenih akcij, ki pa so skupaj kot celota sposobni izvajanja kompleksnih stvaritev (npr. gradnja egipčanskih piramid). Termiti, na primer, so relativno enostavna bitja, ki pa so v interakciji Slika 1. Vzori za razvoj novih vrst algoritmov prihajajo danes iz narave sposobna zgraditi veličastna domovanja in družbe. Najpomembnejši vpliv na razvoj novih vrst algoritmov je (termitnjake). Mravlje in čebele živijo skupaj imela biologija, ki je vplivala na nastanek algoritmov na osnovi zaradi iskanja hrane. Pri tem se mravlje obnašanja inteligence rojev in obnašanja bioloških sistemov. Na razvoj sporazumevajo med seboj posredno s pomočjo novih vrst algoritmov pa vplivajo tudi: kemija, fizika in ostali (družbeni) kemične snovi feromona, katere količina določa vzori. najkrajšo pot do hrane, medtem ko se čebele 3.1. Algoritmi na osnovi inteligence rojev sporazumevajo neposredno prek plesa z mahanjem (angl. waggle dance), s katerim izvidnik (angl. scout) vabi ostale Algoritme na osnovi inteligence rojev odlikujejo naslednje čebele delavke na področja bogata s hrano, t.j. nektarjem. lastnosti: Na razvoj novih vrst algoritmov so vplivali tudi vzori iz - nadzor je porazdeljen med posameznike (decentralizacija), drugih naravoslovnih znanosti, kot npr. kemije in fizike. Kemija - komunikacija poteka lokalno (kolektivno obnašanje), je še posebej vplivala na razvoj t.i. optimizacije na osnovi - vedenje posameznika je podrejeno vedenju sistema (samokemijskih reakcij (angl. Chemical Reaction Optimization, krajše organizacija), CRO), ki izkorišča interakcijo molekul v teh reakcijah za doseganje najnižjega energetskega stanja (Lam & Li, 2012). S - celoten odziv sistema je zelo robusten in prilagodljiv na področja fizike je na razvoj novih vrst algoritmov vplival pojav spremembe v okolju (robustnost in fleksibilnost). simuliranega ohlajanja (angl. Simmulated Annealing, krajše Inteligenca rojev se nanaša na raziskovalno področje, ki se SA), pri katerem tekočino počasi ohlajamo, dokler se snov, v ukvarja s kolektivnim obnašanjem v samo-organiziranih in stanju z minimalno energijo, ne strdi (npr. voda postane led) decentraliziranih sistemih. Ta izraz naj bi prvi uporabil Beni (Kirkpatrick, Gellat, & Veechi, 1983). (1989) pri razvoju celičnih robotov, ki sestojijo iz enostavnih Tudi področje družboslovnih znanosti je dalo svoj pečat pri agentov in se sporazumevajo s pomočjo interakcije z agenti v razvoju novih vrst algoritmov. V tem primeru postane glavni okolici. Danes uporabljamo metode inteligence rojev v vzor človek z različnimi vidiki socialnega vedenja, kot npr. optimizaciji, nadzoru robotov, in usmerjanju v novih procesi zbiranja in soočanja idej (angl. brain-storming), različni generacijah mobilnih omrežij, ki zahtevajo robustnost in čustveni in socialni procesi v družbi in celo procesi v anarhični fleksibilnost. družbi. V glasbi, na primer, je dirigentov cilj poiskati v orkestru Najpomembnejše vrste algoritmov po vzorih iz obnašanja z izvajalci na različnih instrumentih najboljšo harmonijo. Ta inteligentnih rojev so (slika 2): proces je postal vzor za razvoj algoritma iskanja harmonij (angl. Harmony Search, HA) (Geem, 20120). - umetni imunski sistemi (angl. Artificial Immune Systems, krajše AIS) (Dasgupta 1999), 3. Klasifikacija algoritmov po vzorih iz narave - optimizacija z roji delcev (angl. Particle Swarm Optimization, krajše PSO) (Kennedy in Eberhart 1999), Algoritme po vzorih iz narave lahko razdelimo v štiri večje skupine, t.j. algoritme, ki temeljijo na (slika 1): - obnašanju inteligence rojev, naravni evoluciji, fizikalnih in kemijskih sistemih in drugih vzorih iz narave in družbe. V nadaljevanju se osredotočamo na prve tri vrste algoritmov, ki temeljijo na vzorih iz narave. - algoritmi na osnovi cvetnega opraševanja (angl. Flower Pollination, krajše FPA) (Yang, 2012), - algoritmi na osnovi obnašanja netopirjev (angl. Bat Algorithm, krajše BA) (Yang, 2010), - kukavičje iskanje (angl. Cuckoo Search, krajše CS) (Yang in Deb, 2009), - algoritmi na osnovi obnašanja kresnic (angl. Firefly Algorithm, krajše FA) (Yang, 2008), - optimizacija s kolonijami umetnih čebel (angl. Artificial Bee Colony, krajše ABC) (Karaboga & Bastruk, 2007), - optimizacija s kolonijami mravelj (angl. Ant Colony Optimization, krajše ACO) (Dorigo & Di Caro, 1999). Pomurska obzorja 1/ 2014/ 2 | 77 Iztok FISTER et al.: ALGORITMI PO VZORIH IZ NARAVE 3.2. Algoritmi na osnovah Darwinove evolucije Algoritmi temelječi na osnovah Darwinove evolucijske teorije simulirajo delovanje naravne selekcije, ki jo je postavil Darwin (1859). Na podlagi te teorije je nastalo področje evolucijskega računanja (angl. Evolutionary Computation, krajše EC) (Eiben & Smith, 2003). Evolucijsko računanje je sodoben pojem, ki zajema vse algoritme nastale na Darwinovem principu naravne selekcije. Po tem principu imajo za preživetje v naravi največ možnosti najuspešnejši posamezniki, ki prenašajo svoje dobre lastnosti na svoje potomce v procesu reprodukcije (križanje in mutacija). Vse algoritme, ki so nastali na področju evolucijskega računanja, imenujemo tudi evolucijske (angl. Evolutionary Algorithms, krajše EA) in jih delimo v naslednje vrste (slika 3): - genetski algoritmi (angl. Genetic Algorithms, krajše GA) (Goldberg, 1996), Slika 2. Algoritmi na osnovi inteligence rojev spadajo v disciplino umetne inteligence in iščejo vzore za delovanje v kolektivnem obnašanju kolonij mravelj in termitov, rojev čebel in črvov, jatah ptičev in rib. Na sliki so predstavljeni naslednji vzori za nastanek najpogostejših algoritmov, ki si sledijo v smeri urinega kazalca: naravni imunski sistemi, jate ptičev, cvetno opraševanje, obnašanje netopirjev (eholokacija), kukavic (podtakniti svoje jajce v tuje gnezdo), kresnic (svetilnost), kolonij čebel (nabiranje nektarja), in mravelj (feromon). - genetsko programiranje (angl. Genetic Programming, krajše GP) (Koza, 1994), - evolucijske strategije (angl. Evolution Strategies, krajše ES) (Bäck, 1996), - evolucijsko programiranje (angl. Evolutionary Programming, krajše EP) (Fogel, Owens, & Walsh, 1966), - diferencialna evolucija (angl. Differential Evolution, krajše DE) (Storn & Price, 1997) (Brest, Greiner, Bošković, Mernik, & Žumer, 2006). Algoritmi na osnovi inteligence rojev so populacijski, kar pomeni, da namesto ene rešitve vzdržujejo populacijo rešitev, s katerimi simulirajo kolektivno obnašanje delcev (angl. swarm) v roju. Čeprav so vsi algoritmov te vrste zelo podobni med seboj pa se razlikujejo v načinu preiskovanja prostora rešitev, t.j. kako se iz obstoječe pozicije v tem prostoru premakniti na novo, po možnosti Slika 3. Tudi evolucijski algoritmi spadajo v področje umetne boljšo pozicijo. Novo pozicijo v prostoru inteligence. Vsak evolucijski algoritem je sestavljen iz: inicializacije, izračunavamo na podlagi matematičnih enačb, ki selekcije staršev, reprodukcije, ocenitvene funkcije, selekcije preživelih in pogoja ustavljanja. opisujejo pojav, značilen za kolektivno obnašanje vzora, ki ga simuliramo (npr. jate ptičev, kresnic, Čeprav se vse omenjene vrste algoritmov razvijajo ipd.). Običajno se vse rešitve premikajo v smeri trenutne neodvisno druga od druge, pa jih povezujejo podobne lastnosti najboljše in pri tem odkrivajo nove najboljše rešitve. Problem se pri reševanju problemov. Tudi evolucijski algoritmi so v pojavi, ko najboljše rešitve ni mogoče več izboljšati. V tem splošnem populacijski. Vsaka rešitev oz. posameznik v primeru se pojavi stagnacija, ki jo običajno poskušamo reševati populaciji sestoji iz elementov, ki jih v evolucijski terminologiji z dodatnimi prijemi, kot npr. lokalnim iskanjem, ipd… imenujemo tudi gene. V vsakem evolucijskem ciklu so Pri tem je potrebno poudariti, da smo tukaj zajeli kvečjemu posamezniki v populaciji podvrženi procesu reprodukcije, t.j. 10% različnih vrst algoritmov s področja inteligence rojev. delovanju operatorjev križanja in mutacije. Križanje običajno iz Vendar se razvoj novih vrst s tega področja še zdaleč ni končal. dveh posameznikov generira dva potomca, medtem ko mutacija Dan za dnem se pojavljajo nove vrste algoritmov temelječe na spreminja naključno izbran element posameznika. Vsakega inteligenci rojev (npr. družine volkov, psov, mačk, japonskih potomca po reprodukciji ovrednotimo. Pri tem uporabimo žab, ipd…) tako, da se lahko z upravičenostjo sprašujemo, funkcijo uspešnosti povezano s problemom, ki ga rešujemo. V koliko novega prinašajo in ali niso vse te nove vrste algoritmov naslednjo generacijo evolucijskega procesa s selekcijo samo podmnožica že odkritih. Pravo resnico o tem bo prinesel preživelih izberemo najboljše posameznike glede na vrednost čas. funkcije uspešnosti. Evolucijski algoritmi se med seboj razlikujejo glede na predstavitev posameznikov (rešitev). Genetski algoritmi, na primer, delujejo s populacijo posameznikov predstavljenih kot 78 | Pomurska obzorja 1/ 2014/ 2 Iztok FISTER et al.: ALGORITMI PO VZORIH IZ NARAVE dvojiška števila, evolucijske strategije kot realna števila, genetsko programiranje z drevesi implementiranimi v Lispu, in evolucijsko programiranje s končnimi avtomati. Te algoritme danes uspešno uporabljamo na različnih področjih optimizacije, modeliranja in simulacije. 3.3. Algoritmi temelječi na fizikalnih in kemijskih procesih Nekateri današnji algoritmi svojih vzorov ne iščejo v bioloških sistemih, ampak iščejo svoj navdih v fiziki in kemiji. Pri delovanju simulirajo običajno določene fizikalne oz. kemijske zakone, ki vključujejo npr. električni naboj, gravitacijo, rečne sisteme, ipd. Čeprav večina teh vzorov še vedno prihaja iz narave pa jih zaradi svojih posebnosti klasificiramo posebej. Podobno kot algoritmi na osnovi inteligence rojev tudi ti algoritmi iščejo nove rešitve v prostoru preiskovanja po določenih fizikalnih in kemijskih enačbah. Najpomembnejše algoritme te vrste prikazuje slika 4. Celoten spisek algoritmov je zbran v članku Fister Jr., & drugi (2013). Prav tako ti algoritmi niso nujno populacijski. Simulirano ohlajanje (Kirkpatrick, Gellat, & Veechi, 1983), na primer, išče eno samo optimalno rešitev, medtem ko so algoritmi na osnovi inteligentnega odlaganja materiala rek (angl. Intelligent Water Drops, krajše IWD) (Shah-Hosseini, 2007), črnih lukenj (angl. Black Holes, krajše BH) (Hatamlou, 2012), galaktičnega iskanja (angl. Galaxy-based Search Algorithm, krajše GbSA) (ShahHoseini, 2013), kemijskih reakcij (Lam & Li, 2012) in dogajanja v obremenjenih sistemih (angl. Charged System Search, krajše CSS) (Kaveh & Talatahari, 2010) populacijski. Med to vrsto algoritmov lahko najdemo tudi večagentne sisteme, kot npr. gravitacijski algoritem (angl. Gravitational Search Algorithm, krajše GSA) (Rashedi, Nezamabadi-pour, & Saryazdi, 2009). 4. Povzetek Cilj tega podpoglavja je narediti kratek povzetek klasifikacije algoritmov po vzoru iz narave, omeniti njihove avtorje in dati reference na njihova dela. S tem želimo pomagati raziskovalcem oz. razvijalcem programske opreme, po eni strani, pri iskanju orodij za reševanje njihovih problemov in jih, po drugi strani, spodbuditi za uporabo algoritmov po vzoru iz narave. Vzore iz narave, pripadajoče algoritme, njihove oznake in avtorje prikazuje tabela 1. Slika 4. Osnova algoritma SA predstavlja enačba ohlajanja, ki spreminja verjetnost premikanja iz ene točke v prostoru preiskovanja v drugo. Algoritem IWD izkorišča odlaganje materiala v rekah, ki vpliva na optimalni tok reke brez nepotrebnih okljuk (meandrov). BH temelji na pojavu črnih lukenj, ki absorbirajo najslabše posameznike v populaciji. GbSA deluje na osnovi pojava spiralnih rok pri nekaterih galaksijah, pri čemer išče njihove mejne vrednosti. GSA temelji na Newtonovem zakonu o gravitaciji. CRO pri optimizaciji simulira delovanje prvih dveh zakonov termodinamike, ki spremljata kemijske reakcije, medtem ko CSS izkorišča Coulombov zakon elektrostatike in Newtonove zakone mehanike. Pomurska obzorja 1/ 2014/ 2 | 79 Iztok FISTER et al.: ALGORITMI PO VZORIH IZ NARAVE Tabela 1. Algoritmi iz narave. Vzor iz narave Algoritem Oznaka Avtor Imunski sistemi Umetni imunski sistemi AIS Dasgupta (1999) Roji delcev Optimizacija z roji delcev PSO Kennedy in Eberhart (1999) Cvetno opraševanje Alg. na osnovi cvetnega opraševanja FPA Yang (2012) Eholokacija netopirjev Alg. na osnovi obnašanja netopirjev BA Yang (2010) Podtikanje jajc Kukavičje iskanje CS Yang in Deb (2009) Svetilnost Alg. na osnovi obnašanja kresnic FA Yang, 2008 Iskanje hrane Optimizacija z kolonijami umetnih čebel ABC Karaboga & Bastruk (2007) Darwinova evolucija Genetski algoritmi GA Goldberg (1996) Darwinova evolucija Genetsko programiranje GP Koza (1994) Darwinova evolucija Evolucijske strategije ES Bäck (1996) Darwinova evolucija Evolucijsko programiranje EP Fogel, Owens, & Walsh (1966) Darwinova evolucija Diferencialna evolucija DE Storn & Price (1997) Simulirano ohlajanje Alg. na osnovi simuliranega ohlajanja SA Kirkpatrick, Gellat, & Veechi (1983) Odlaganje prsti rek Inteligentno odlaganje materiala rek IWD Shah-Hosseini (2007) Črne luknje Alg. na osnovi črnih lukenj BH Hatamlou (2012) Galaksije Galaktično iskanje GbSA Shah-Hoseini (2013) Kemijske reakcije Optimizacija s kemijskimi reakcijami CRO Lam & Li (2012) Newtonov zakon Gravitacijski iskalni algoritem GSA Rashedi, Nezamabadi-pour, & Saryazdi (2009) gravitacijski Ime algoritmi po vzoru iz narave ni nastalo naključno. V Sloveniji smo že leta 2004 ustanovili skupino z imenom Algoritmi po vzoru iz narave (AVN), ki organizira delavnice s področja evolucijskega računanja in drugih računalniških metod po vzorih iz narave. Na delavnicah se srečujejo podiplomski študenti in mentorji iz raziskovalnih ustanov, univerz in podjetij. Delavnice potekajo izmenično na Institutu Jožef Stefan v Ljubljani in Fakulteti za elektrotehniko, računalništvo in informatiko v Mariboru. K udeležbi so vabljeni vsi, ki se ukvarjajo s temi metodami in njihovo uporabo in želijo aktivno sodelovati s svojimi prispevki ter pridobivati in deliti izkušnje na tem področju. kemijskih zakonih. V ta pregled nismo zajeli algoritmov, ki simulirajo delovanje človeka in njegovih procesov pri interakciji z drugimi ljudmi. V nadaljevanju dela se želimo posvetiti tudi temu področju. Viri 1. Bäck, T. (1996). Evolutionary Algorithms in Theory and Practice Evolution Strategies, Evolutionary Programming, Genetic Algorithms. Oxford University Press. 2. Beni, G., & Wang, J. (1989). Swarm Intelligence in Cellular Robotic Systems. Proceedings of NATO Advanced Workshop on Robots and Biological Systems, 26-30. Tuscany, Italy. 3. Blum, C., & Merkle, D. (2008). Swarm Intelligence. Berlin: Springer-Verlag. 4. Brest, J., Greiner, S., Bošković, B., Mernik, M., & Žumer, V. (2006). Self-adapting Control Parameters in Differential Evolution: A comparative Study on Numerical Benchmark Problems. IEEE Transaction on Evolutionary Computation , 10(6), 646-657. 5. Darwin, C. (1859). On the Origin of Species. London: Harvard University Press (1964). 5. Zaključek V tem članku predstavljamo novo klasifikacijo algoritmov po vzoru iz narave, da bi pomagali tistim, ki s temi metodami želijo reševati vsakodnevne probleme v praksi. To je področje umetne inteligence, ki je dandanes izredno dinamično. Skoraj vsak dan se namreč pojavljajo nove vrste algoritmov, zato njihova klasifikacija še ni določena. Po eni strani jih računalniški strokovnjaki uvrščajo med evolucijske algoritme, po drugi strani med mehko računanje (angl. soft computing). Sami zastopamo stališče, da gre v tem primeru za algoritme, ki jih lahko uvrščamo v t.i. računsko inteligenco (angl. Computational Intelligence), katera išče vzore za svoje delovanje v naravi, kot npr. v Darwinovem zakonu boja za obstanek, obnašanju različnih vrst žuželk, živali in rastlin, ter v fizikalnih in 80 | Pomurska obzorja 1/ 2014/ 2 Iztok FISTER et al.: ALGORITMI PO VZORIH IZ NARAVE 6. Das, S., & Suganthan, P. (2011). Differential Evolution: A Survey of the State-of-the-art. IEEE Transaction on Evolutionary Computation , 15(1), 4-31. 21. Koza, J. (1994). Genetic Programming 2 - Automatic Discovery of Reusable Programs. Cambridge, USA: MIT Press. 7. Dasgupta, D. (1999). Information Processing in the Immune System. V D. Corne, M. Dorigo, & F. Glover, New Ideas in Optimization, 161-167. New York, USA: McGrawHill. 22. Lam, A. Y., & Li, V. O. (2012). Chemical Reaction Optimization: A Tutorial. Memetic Computing , 3-17. 8. 9. Dorigo, M., & Di Caro, G. (1999). The Ant Colony Optimization Meta-heuristic. V D. Corne, M. Dorigo, & F. Glover, New Ideas in Optimization, 11-32. London, UK: McGraw Hill. Eiben, A., & Smith, J. (2003). Introduction to Evolutionary Computing. Berlin: Springer-Verlag. 10. Fister, I. Jr., Fister, I., Brest, J., & Fister, D. (2013). A Brief Review of Nature-inspired Algorithms for Optimization. Electrotechnical Review , 80(3), 161-168. 11. Fister, I., Fister, I. Jr., Brest, J., & Žumer, V. (2012). Memetic Artificial Bee Colony Algorithm for Large-scale Global Optimization. IEEE Congress on Evolutionary Computation, 1-8. 12. Fogel, L., Owens, A., & Walsh, M. (1966). Artificial Intelligence through Simulated Evolution. New York, US: John Willey. 13. Garey, M., & Johnson, D. (1979). Computers and Intractability: A Guide to the Theory of NP-Completeness. New York, NY, USA: W.H. Freeman & Co. 23. Storn, R., & Price, K. (1997). Differential Evolution: A Simple and Efficient Heuristic for Global Optimization over Continuous Spaces. Journal of Global Optimization , 11(4), 341-359. 24. Yang, X.-S. (2010). A New Metaheuristic Bat-inspired Algorithm. V C. Cruz, J. Gonzales, G. Krasnogor, & D. Pelta, Nature Inspired Cooperative Strategies for Optimization (NISCO 2010), 65-74. Springer-Verlag: Berlin. 25. Yang, X.-S. (2008). Firefly Algorithm. V X.-S. Yang, Nature-Inspired Metaheuristic Algorithms, 79-90. London, UK: Luniver Press. 26. Yang, X.-S. (2012). Flower Pollination Algorithm for Global Optimization. V J. Durand-Lose, & N. Jonoska, Unconventional Computation and Natural Computation, 240-249. Berlin: Springer-Verlag. 27. Yang, X.-S., & Deb, S. (2009). Cuckoo Search via Levy Flights. World Congress & Biologically Inspired Computing (NaBIC 2009), 210-214. IEEE Publication. 14. Geem, Z. W. (2012). Recent Advances in Harmony Search Algorithm. Berlin: Springer-Verlag. 15. Goldberg, D. (1996). Genetic Algorithms in Search, Optimization, and Machine Learning. MA: AddisonWesley. 16. Holland, J. (1992). Adaptation in Natural and Artificial Intelligence: An Introductory Analysis with Applications to Biology, Control and Artificial Intelligence. Cambridge, MA, USA: MIT Press. 17. Karaboga, D., & Bastruk, B. (2007). A Powerful and Efficient Algorithm for Numerical Function Optimization: Artificial Bee Colony (ABC) Algorithm. Journal of Global Optimization , 39(3), 459-471. 18. Kennedy, J., & Eberhart, R. (1999). The Particle Swarm Optimization; Social Adaptation in Information Processing. V D. Corne, M. Dorigo, & F. Glover, New Ideas in Optimization, 379-387. London, UK: McGraw Hill. 19. Kirkpatrick, S., Gellat, C. J., & Veechi, M. (1983). Optimization by Simulated Annealing. Science , 220(4578), 671-680. 20. Korošec, P., Šilc, J., & Filipič, B. (2012). The Differential Ant-stigmegy algorithm. Information Sciences , 82-97. Pomurska obzorja 1/ 2014/ 2 | 81 Tehnika Matej Zadravec1, Uroš Jeke2, Matjaž Hriberšek1 Numerična simulacija uparjanja kapljevitega filma POVZETEK V večini tehniških in naravnih sistemov je prisoten vlažen zrak. Vlažen zrak lahko pri določenem tlaku in temperaturi sprejme določeno količino vlage, pri čemer pride do točke nasičenja. Kondenzacija na stenah trdnih površin velikokrat ovira normalno delovanje v tehničnih napravah in je nezaželena. Nasproten pojav kondenzacije je uparjanje. Predstavljen je splošen matematičnofizikalni model uparjanja kapljevitega filma in njegova implementacija v obstoječi programski paket za računalniško dinamiko tekočin. Implementacija numeričnega modela uparjanja je bila preverjena na primeru simulacije pravokotnega 3D kanala in primerjana z izvedenim eksperimentom uparjanja kapljevitega filma. Zrak, ki teče nad vodno površino uparja vodo, pri čemer se vodna površina ohlaja, dokler se ne vzpostavi ravnotežje med prevodom toplote iz zraka in latentne toplote, ki se porabi za uparjanje. Spremljali smo kombiniran prenos toplote in snovi. Ugotovljeno je bilo, da je vpliv naravne konvekcije v primerjavi s konvektivnim tokom zanemarljiv in da je bila implementacija modela uparjanja uspešna. Ključne besede: uparjane, računalniška dinamika tekočin, naravna konvekcija, mešana konvekcija. 1. Uvod V mnogih inženirskih problemih predvsem pri procesih klimatizacije in prezračevanja in marsikaterem procesnem postrojenju se srečujemo s procesom uparjanja in kondenzacije. V vseh teh inženirskih problemih imamo večinoma prisoten vlažen zrak, ki pa vemo vsebuje pri nekem stanju temperature in tlaka določeno količino vlage. Tak vlažen zrak je sposoben sprejeti še dodatno količino vlage vse do točke nasičenja vlažnega zraka, pri kateri se začne vlaga izločati iz vlažnega zraka v obliki kapljevite faze in sicer v obliki kapljic. Točka nasičenja je odvisna od temperature in tlaka pri kateri se nahaja vlažen zrak in zato zasledimo plast kondenzirane vlage na hladnejših površinah, kjer je na teh površinah njihova temperatura nižja od temperature rosišča pri določenem tlaku. Vsa ta stanja je zelo dobro možno opisati s pomočjo analitičnih in empiričnih enačb za enostavne primere, ko imamo znane hitrosti in ostale fizikalne parametre sistema, pri čemer pa je uporaba teh analitično-empiričnih nastavkov v bolj geometrijsko in fizikalno kompleksnejših geometrijah nemogoča zaradi različnih vplivov na sam proces uparjanja in kondenzacije. Dandanes obstaja precej računalniških paketov, ki so sposobni reševati kompleksne inženirske problema s področja dinamike tekočin. Ampak vseh fizikalnih pojavov seveda nimajo vključenih in zajetih, kar pa je seveda možno v same programske pakete vgraditi preko zunanjih računalniških podprogramov, ki jih lahko uporabnik sprogramira in implementira v obstoječi računski paket. V našem primeru reševanja prej omenjenih fizikalnih pojavov smo se odločili za začetni korak implementacije procesa uparjanja preko različnih empiričnih nastavkov izračuna stanj in točke nasičenja vlažnega zraka. Kot testni primer je bil izbran preprost primer izhlapevanja vode iz plasti kapljevine v kanalu, skozi katerega teče zrak. Primer je dobro obdelan eksprimentalno kot tudi numerično v članku Talukdarja [1]. Predstavljen prispevek se navezuje na na magistrsko delo Jekeja [2]. 2. Definicija problema Vključeni teoretični model uparjanja je preizkušen na modelu kratkega 3D pravokotnega kanala prikazanega na sliki 1, katerega geometrija je bila povzeta po izhodiščnem delu [1]. Simulacija je bila izvedena za vrednost Reynoldsovega števila Re=844. Slika 1. – Geometrija numeričnega modela [2]. 1 82 | Fakulteta za strojništvo, Univerza v Mariboru 2 Hella Saturnus Slovenija d.o.o. Pomurska obzorja 1/ 2014/ 2 Matej ZADRAVEC, Uroš JEKE, Matjaž HRIBERŠEK: NUMERIČNA SIMULACIJA… Zaradi zahtevnosti obravnavane problematike bomo modelirali enofazni dvokomponentni tok, pri čemer se bo druga komponenta zaradi majhnih masnih deležev obravnavala ločeno od zraka in sicer s posebno enačbo ohranitve vodne pare v obliki koncentracije vodne pare. Sistem vodilnih enačb, ki jih moramo rešiti, da dobimo podatke o lokalnem koncentracijskem polju vodne pare kot tudi o difuzijskih tokovih na robu območja, je razširjeni sistem Navier-stokes enačb (en.1-6): 𝜕𝑢 𝜕𝑥 [𝑢 + + 𝜕𝑤 𝜕𝑧 = 0, (1) 𝜕𝑢 𝜕𝑢 𝜕𝑢 +𝑣 +𝑤 ]= 𝜕𝑥 𝜕𝑦 𝜕𝑧 = − [𝑢 𝜕𝑣 𝜕𝑦 1 𝜕𝑝 𝜌0 𝜕𝑥 + 𝜈0 [ 𝜕2𝑢 𝜕𝑥 2 + 𝜕2 𝑢 𝜕𝑦 2 (2) + 𝜕2𝑢 𝜕𝑧 2 ], 𝜕𝑣 𝜕𝑣 𝜕𝑣 1 𝜕𝑝 𝜕 2 𝑣 𝜕 2𝑣 𝜕 2𝑣 +𝑣 +𝑤 ]= − + 𝜈0 [ 2 + 2 + 2 ] + 𝜕𝑥 𝜕𝑦 𝜕𝑧 𝜌0 𝜕𝑦 𝜕𝑥 𝜕𝑦 𝜕𝑧 (3) + 𝛽𝑇 𝑔(𝑇 − 𝑇0 ) + 𝛽𝐶 𝑔(𝐶 − 𝐶0 ), [𝑢 𝜕𝑤 𝜕𝑥 +𝑣 𝜕𝑇 𝜕𝑤 𝜕𝑦 +𝑤 𝜕𝑇 𝜕𝑤 𝜕𝑧 ]= − 1 𝜕𝑝 𝜌0 𝜕𝑧 + 𝜈0 [ 𝜕 2𝑇 𝜕𝑇 𝜕2 𝑤 𝜕𝑥 2 𝜕 2𝑇 + 𝜕2𝑤 𝜕𝑦 2 + 𝜕2𝑤 𝜕𝑧 2 𝜕 2𝑇 [𝑢 𝜕𝑥 + 𝑣 𝜕𝑦 + 𝑤 𝜕𝑧 ] = 𝜆0 [𝜕𝑥 2 + 𝜕𝑦2 + 𝜕𝑧 2 ], 𝜕𝐶 𝜕𝐶 𝜕 2𝐶 𝜕𝐶 𝜕 2𝐶 𝜕 2𝐶 [𝑢 𝜕𝑥 + 𝑣 𝜕𝑦 + 𝑤 𝜕𝑧 ] = 𝐷𝑣 [𝜕𝑥 2 + 𝜕𝑦2 + 𝜕𝑧 2 ] + 𝑟. ], 1 𝜕𝑝 1 1 0 𝛽𝐶 = − = 1 𝜕𝑝 1 𝑀𝑣 = [ − 1] = 𝜌 𝜕𝐶 𝜌0 𝑀𝑧 1∙𝑚3 [ 28,960 𝑘𝑔/𝑘𝑚𝑜𝑙 1,185 𝑘𝑔 18,016 𝑘𝑔/𝑘𝑚𝑜𝑙 𝑁𝑢𝑧 = 𝛼∙𝐷ℎ 𝑆ℎ𝑧 = 𝐾∙𝐷ℎ − = 𝜆 = 𝐷 𝜕𝑇 | 𝜕𝑦 𝑦=−0,01025 𝑇𝑤 − 𝑇𝑚 − 𝜕𝐶 | 𝜕𝑦 𝑦=−0,01025 𝐶𝑤 − 𝐶𝑚 𝐷ℎ , (9) 𝐷ℎ . (10) kjer se indeks w nanaša na temperaturo oz. koncentracijo na vodni površini in indeks m na povprečno temperaturo oz. koncentracijo po celotnem volumnu kanala. Položaj y=-0,01025 je izbran tako, da sovpada s položajem površine vodne gladine. Na površini vodne gladine je za izračun parcialnega tlaka nasičenja vodne pare uporabljena Antoine-jeva enačba (en.11), ki velja za tlak 𝑝 = 1 𝑏𝑎𝑟 in temperaturno območje 274 𝐾 < 𝑇 < 373 𝐾: (4) 3984,923 𝑝𝑣,𝑠 = 𝑒𝑥𝑝 (11,96481 − 𝑇−39,724) [𝑏𝑎𝑟], (11) (5) (6) Upoštevali bomo zgolj časovno ustaljeno stanje, seveda pa je model možno razširiti za simulacijo časovno spremenljivih razmer. Vpliv prisotnosti vodne pare v zraku je upoštevan z dodatnimi vzgonskimi silami v enačbi ohranitve gibalne količine (en.3), temperaturni razteznostni koeficient 𝛽𝑇 in snovni razteznostni koeficient 𝛽𝐶 podamo z sledečima enačbama (en.7 in en.8) 𝛽𝑇 = − 𝜌 𝜕𝑇 = 𝑇 = 298,15 𝐾 = 0,00335 𝐾 −1, Upoštevan bo tudi prenos toplote med zrakom in kapljevinskim filmom. Na vodni gladini bomo spremljali snovni uparjalni tok, ter lokalne in povprečne vrednosti Nusselt-ovega ter Sherwoodovega brezdimenzijskega števila (en.9 in en.10): (7) (8) − 1] = 0,513 𝑚 3 /𝑘𝑔, kjer sta 𝑀𝑧 in 𝑀𝑣 molski masi zraka in vode. Pri termičnem razteznostnem koeficientu upoštevamo temperaturo, pri kateri bodo upoštevane lastnosti zraka v numeričnem modelu, in sicer 25°𝐶, kar pomeni, da je bila vrednost gostote 𝜌0 = 1,185 𝑘𝑔/ 𝑚 3. Za sistem zrak-voda smo uporabili vrednost snovne difuzivnosti 𝐷𝑣 = 2,82 ∙ 10−5 𝑚2 /𝑠 [3]. Snovne lastnosti zraka so bile konstantne. Uporabljene lastnosti zraka so bile pri temperaturi 𝑇 = 25 °𝐶 in referenčnem tlaku 𝑝 = 1 𝑏𝑎𝑟. Vpliv vlage na gibanje zraka bo upoštevan s pomočjo dodatne konvektivne hitrosti na meji med vodno površino in zrakom, ki je posledica procesa izhlapevanja vode. Predpostavili bomo, da je glavni mehanizem prenosa snovi v tanki mejni plasti tik ob vodni površini difuzija, zato lahko upoštevamo model enostranske difuzije – Stefanov tok, [4]. Uparjanje bomo predpisali kot tok z neko hitrostjo glede na normalo na površino. Koncentracijo vodne pare v zraku (en.12) naposled izračunamo s pomočjo definicije masnega deleža (en.13) in definicije vlažnosti zraka (en.14) 𝐶= 𝑚𝑣 𝜉𝑣 = 𝑉 𝑚𝑣 𝑚 (12) = 𝜉𝑣 𝜌, =𝑚 𝑥 = 0,622 ∙ 𝑚𝑣 𝑣 +𝑚𝑧 𝑝𝑣,𝑠 𝑝−𝑝𝑣,𝑠 𝑥 = 1+𝑥, (13) = 𝑥𝑠 (𝑝, 𝑇). (14) 3. Numerični model Numerični model kanala smo diskretizirali s strukturirano heksaedersko mrežo, ki je imela 105.896 vozlišč kar nam je omogočalo dovolj natančne rezultate ob sprejemljivih računskih časih. Največji izziv pri implementaciji modela uparjanja je bilo direktno predpisovanje robnih pogojev odvoda spremenljivke, bodisi po komponenti kartezijevega koordinatnega sistema bodisi po normali na površino. Rešitev upoštevanja takih robnih pogojev na vodni površini smo našli v implementaciji podprograma napisanega v programskem jeziku Fortran 90, ki se izračuna izven osnovnega računskega modula reševanja sistemov enačb in se potem podatki preko vmesnika prenašajo iz osnovnega računskega algoritma v ta podprogram. Osnovni koncept je prikazan na sliki (2). Pomurska obzorja 1/ 2014/ 2 | 83 Matej ZADRAVEC, Uroš JEKE, Matjaž HRIBERŠEK: NUMERIČNA SIMULACIJA… ki je enaka hitrosti uparjanja Temperatura na vodni površini bo nižja kot pa je temperatura toka zraka. Namesto konstantnih vrednosti temperature in koncentracij na vodni površini se bosta temperatura in koncentracija izračunavali na podlagi robnega pogoja za energijo, ki pravi, da je latentna toplota uparjanja enaka prenosu toplote med zrakom in vodo. Uparjanje vlage je v tem modelu upoštevano kot tok s hitrostjo 𝑣𝑤 . Hitrost je glede na normalo na površino enaka 𝑣𝑤 = ( 𝜕 𝜉𝑣 𝜕𝑛 steno, torej je vrednost toplotnega toka enaka )) . 𝑠 𝑄̇𝑣 𝐴 = 0. Ob tej predpostavki je na našem robnem pogoju količina toplote potrebna za fazno spremembo pri uparjanju vode enaka senzibilnemu prenosu toplote iz zraka na vodno površino. Robni pogoj za zrak na vodni površini lahko zapišemo kot Robni pogoji ki smo jih predpisali na omenjenem modelu se lahko delijo na štiri osnovne sklope:  Za vtok zraka je bil predpisan analitični hitrostni profil zraka v pravokotnem kanalu [6]. Povprečna hitrost 𝑤𝑎𝑣 je definirana preko Reynoldsovega števila na vstopnem robu, kjer smo upoštevali 𝑅𝑒 = 844, kar nam okarakterizira laminaren tok vlažnega zraka [3]. Temperatura na vstopu je bila izbrana, glede na eksperimentalne in numerične rezultate v izhodiščnem delu [1] in znaša 𝑇0 = 22,4 °𝐶. Na vstopnem robu predpišemo še vrednost koncentracije vodne pare. Ta se izračunava kot t.i. dodatna spremenljivka. Vrednost koncentracije vlage v vlažnem zraku je odvisna od temperature na vstopu 𝑇0 in relativne vlažnosti zraka 𝜑0 = 53,1 % [1].  Vodna površina je definirana kot stena na kateri imamo izvor toplote in dodatno hitrost zraka v domeno s hitrostjo, ( Pri pojavu uparjanja in kondenzacije imamo zaradi spremembe agregatnega stanja intenziven prenos toplote. Ta toplotni tok vpliva na ravnotežno stanje na vodni gladini. Toplotni tok, ki priteka na vodno površino je vsota prevoda toplote v tekočini tik ob steni in latentne toplote, ki se pri uparjanju porabi, pri kondenzaciji pa sprošča. Predpostavimo, da imamo na vodni površini adiabatno Slika 2. Numerični algoritem reševanja s pomočjo zunanje rutine v programski paket Ansys CFX 13.0  Na stranskih stenah in zgornji steni kanala smo predpisali brez zdrsni robni pogoj, kar pomeni da so vrednosti hitrosti zraka enake nič. Na teh stenah je bil predpisan adiabatni robni pogoj. Preko sten ni nobenega toka dodatne spremenljivke za vlago in je odvod koncentracije po normali enak nič. 𝐷𝑣 1−𝜉𝑣 𝜕𝑇 −𝜆 𝜕𝑛 = 𝑟0 ∙ 𝐷𝑣 ∙ 𝜕𝐶 𝜕𝑛 . Vrednost koncentracije vodne pare tik ob vodni površini temelji na temperaturi vode 𝑇𝑤 . Temperatura se izračunava 𝜕𝑇 iz enačbe −𝜆 𝜕𝑛 = 𝑟0 ∙ 𝐷𝑣 ∙ 𝜕𝐶 𝜕𝑛 in predpostavlja razmere popolnoma nasičenega zraka. Vrednost koncentracije vlage 𝑝 izračunamo iz enačbe 𝐶 = 0,622 ∙ 𝑣,𝑠 ∙ 𝜌0 , kjer dobimo 𝑝−𝑝𝑣,𝑠 vrednost 𝑝𝑣,𝑠 iz Antoine-ove enačbe (en.11). Na spodnji površini ne smemo pozabiti, da je širina vodne površine ožja kot je širina kanala. Na preostankih spodnjega roba, ki ni del vodne površine, predpostavljamo adiabaten robni pogoj, kjer so stene nepropustne za toplotni tok in snovni tok vodne pare.  Na izstopnem robu nimamo nobenim gradientov, ki bi bili gonilna sila in predpišemo 𝜕𝑤 𝜕𝑧 𝜕𝑇 sprememb 𝜕𝑢 𝜕𝑥 = 𝜕𝑣 𝜕𝑦 = 𝜕𝐶 = 𝜕𝑧 = 𝜕𝑧 = 0. 4. Rezultati Na vstopnem robu imamo predpisani hitrostni profil kar lahko vidimo na sliki 3 in razvidno je tudi kako se tok zraka v kanalu še razvija in se v 0,6 𝑚 dolžine popolnoma še ne razvije. Slika 3. Hitrost v prerezni ravnini YZ pri x=0,00 m in detajl hitrosti na vstopnem robu. 84 | Pomurska obzorja 1/ 2014/ 2 Matej ZADRAVEC, Uroš JEKE, Matjaž HRIBERŠEK: NUMERIČNA SIMULACIJA… Na vodni površini smo upoštevali mehanizem enostranske difuzije in s tem gibanje vodne pare (uparjanje vode) opisali kot tok s hitrostjo podano na sliki 4 (levo). Vrednosti hitrosti zaradi uparjanja vidimo da so zelo nizke v primerjavi s konvektivno hitrostjo, ki smo jo predpisali toku zraka na vstopu. Za potrditev rezultatov smo izrisali še diagram povprečnih vrednosti hitrosti na vodni površini v prečni smeri na tok (smer x) po dolžini kanala v primeru z upoštevanjem difuzijske hitrosti na medfazni meji in brez upoštevanja. (slika 4 - desno). Če opazujemo krivuljo difuzijske hitrosti ima ta eksponentno obliko. Proti koncu kanala se vrednost približuje konstantni vrednosti, vendar se še ne ustali. To lahko povežemo z razvojem mejne plasti, ki se sodeč po difuzijski hitrosti popolnoma še ne razvije. Vrednosti difuzijske hitrosti so v primerjavi z ostalimi hitrostmi v kanalu nekaj redov manjše, zato je pričakovan zanemarljiv vpliv. Na sliki 5 so prikazane konture temperatur in koncentracijske porazdelitve (slika 6). Iz slik je vidno, da se termična in koncentracijska mejna plast še razvijata. Iz literature je moč zaslediti, da se termična mejna plast popolnoma razvije pri toku z 𝐷ℎ = 0,0384 𝑚, 𝑅𝑒 = 844 in 𝑃𝑟 = 0,71 šele po več kot 1 𝑚. Iz dobljenih rezultatov lahko ugotovimo, da se kljub nižji vrednosti karakterističnega števila koncentracijski profil še ne razvije popolnoma po dolžini kanala 𝐿 = 0,6 𝑚. Rezultati Nusseltovega in Sherwoodovega števila na sliki 7 so podani za primer brez upoštevanja difuzijske hitrosti na vodni površini in z upoštevanjem difuzijske hitrosti ter so primerjani z vrednostmi, ki jih je v svoji simulaciji pridobil Talukdar [1]. Vrednosti so povprečne vrednosti prečno na smer toka (v 𝑥 smeri) na ravnini 𝑋𝑍, kjer je 𝑦 = −0,01025 𝑚. Iz teh rezultatov lahko ugotovimo, da upoštevanje difuzijske hitrosti na vodni površini ne vpliva na končni rezultat. Dobljene vrednosti na vstopnem delu kanala odstopajo od Talukdarjevih rezultatov [1], ko se pomikamo vzdolž kanala (po 𝑧 smeri) se rezultati približujejo in vse bolj ujemajo. Vzrok v delnem razhajanju je v primerjalno redki računski mreži, katere vpliv je največji prav na področju najvišjih snovnih tokov, to je na vstopu v kanal. Slika 4. Hitrosti na vodni površini. Slika 5. Izoterme na ravnini YZ pri x=0,00 m. Slika 6. Konture porazdelitve koncentracij vodne pare na ravnini YZ pri x=0,00 m. Pomurska obzorja 1/ 2014/ 2 | 85 Matej ZADRAVEC, Uroš JEKE, Matjaž HRIBERŠEK: NUMERIČNA SIMULACIJA… NuZ Talukdar [13] ShZ Talukdar [13] NuZ brez dodane hitrosti ShZ brez dodane hitrosti NuZ dodana difuzijska hitrost Nuz, Shz 20,000 18,000 16,000 14,000 12,000 10,000 8,000 6,000 4,000 2,000 ,000 00,001 00,010 z/(Re Pr Dh) Slika 7. Razporeditev lokalnih vrednosti Nusseltovega in Sherwoodovega števila vzdolž kanala. 5. Zaključek Izdelani računalniški podprogram, ki je bil uporabljen v programskem paketu Ansys-CFX [5], je primeren za modeliranje procesa uparjanja na vodnih površinah. Rezultati izvedenih numeričnih simulacij za izbrani testni primer dokazujejo, da se rezultati dobro ujemajo z referenčnimi rezultati [1]. Pri tem je potrebno poudariti, da je za uporabo modela na gostejših računskih mrežah, kjer je vpliv izmenjave toplote tik ob vodni površini lokalno zelo velik, potrebno podrobneje preučiti način vključitve te vrste prenosa toplote v računski algoritem Ansys-CFX, predvsem s ciljem, da zagotovimo monotono konvergenco celotnega sklopa enačb. Literatura 1. Talukdar Prabal, Iskra Conrad R., Simonson Carey J. Combined heat and mass transfer for laminar flow of moist air in a 3D rectangular duct: CFD simulation and validation with experimental data. International Journal of Heat and Mass Transfer, 2008 , 51 , 3091 – 3102. 2. U. Jeke. Modeliranje kondenzacije in uparjanja vlage na trdnih površinah z računalniško dinamiko tekočin. Magistrsko delo študijskega programa 2. stopnje Strojništvo. Fakulteta za strojništvo, Univerza v Mariboru. Avgust 2012. 3. Škerget Leopold. Prenosni pojavi: Prenos gibalne količine, toplote in snovi. Osnutek učbenika., 2012. 4. L. Škerget, M. Hriberšek, J. Ravnik, M. Ramšak. Modeliranje kondenzacije in uparjanja vlage na trdnih površinah svetil z računalniško dinamiko tekočin. Poročilo, Univerza v Mariboru, Fakulteta za strojništvo, Inštitut za energetsko, procesno in okoljsko inženirstvo, 2010. 5. Ansys inc. Ansys CFX 13.0, 2010. 6. Shah R.K., London A. Chapter VII; Rectangular Ducts. V: Shah R.K., London A. Advances in Heat Transfer: L. Laminar Flow Forced Convection in Ducts. New York: Academic Press, 1978. 86 | Pomurska obzorja 1/ 2014/ 2 00,100 Tehnika Matej Zadravec1,*, Matjaž Hriberšek1, Jure Ravnik1 Ločevanje magnetnih delcev iz tekočine POVZETEK V številnih vejah znanosti je danes možno zaslediti uporabo nano in mikro delcev. V prispevku so obravnavani magnetni delci v tekočini pod vplivom sile neuniformnega magnetnega polja. Za obravnavo omenjenega problema je bila v sklopu računalniške dinamike tekočin uporabljena metoda robnih elementov. Uporabljene enačbe temeljijo na hitrostno vrtinčni formulaciji NavierStokesovih enačb za izračun toka tekočine in Langrangeovega pristopa pri sledenju delcev. Obravnavana je enosmerna interakcija med tekočino in delci. Izpeljani numerični algoritem se uporabi pri študiju vpliva geometrije izločevalnega kanala, predvsem v notranjosti kanala nameščene pregrade, na samo učinkovitost izločevanja delcev. Umeščene trapezoidne pregrade v kanalu pravokotne oblike oblikujejo bifurkacije v samem kanalu, kar bi naj povečalo učinkovitost izločanja magnetnih delcev na stene kanala zaradi vpliva magnetnih in hidrodinamskih sil. Učinkovitost izločanja je bila določena za različne jakosti magnetnega polja in pretoke tekočine. Poleg analize vpliva pretočnih vrednosti in magnetnih sil je bila narejena tudi analiza vpliva velikosti pregrad. Rezultati so se primerjali z rezultati izločanja v ravnem kanalu brez pregrad. Rezultati analize kažejo na smiselnost uporabe manjših pregrad v primeru nizkih pretočnih hitrosti. V primeru višjih pretočnih hitrosti pa se je izkazalo izločanje delcev kot manj učinkovito v primerjavi z izločanjem v ravnem kanalu brez pregrad. Ključne besede: sledenje delcev, računalniška dinamika tekočin, metoda robnih elementov, ločevanje pod vplivom magnetnega polja. 1. Uvod V naravi zasledimo mnogo sistemov v katerih se pojavljajo suspenzije. Suspenzije so sestavljene iz osnovnega in suspendiranega medija, v večini primerov so to zmesi tekočine in trdnih delcev. Ločevalne naprave imajo nalogo iz suspenzije odstraniti trdne delce. Te naprave delujejo na principu izkoriščanja fizikalnih lastnosti delcev in jih ločujejo npr. na podlagi gravitacijskih sil (usedalniki), centrifugalnih sil (cikloni), in v primeru mikro in nano delcev na podlagi magnetnih sil. Ločevalne naprave s pomočjo magnetnega polja (High Gradient Magnetic Separation devices-HGMS), ki so bili uporabljene v dosedanjih raziskavah [1-5], delujejo na principu ustvarjenega magnetnega polja s pomočjo dveh vodnikov. Uporaba takšnih sistemov je prisotna na številnih področjih znanosti v naravoslovju, medicini in tehniških ved. V biomedicini se na primer uporabljajo za magnetno separacijo želenih celic in drugih bioloških enot [6-8] . Poleg tega se uporabljajo delci z magnetnimi lastnostmi, za dostavo terapevtskih drog na določena mesta v telesu [9]. Številni raziskovalci [10-12] so usmerili tudi svoje raziskave v sintezo delcev z ugodnimi magnetnimi lastnostmi in biološko združljivostjo. V tem članku bo opravljena numerična simulacija ločevanja toka z vključenimi magnetnimi delci pod vplivom zunanjega magnetnega polja na primeru ravnega kanala z vstavljeno 1 pregrado. Uporabljena bo 3D numerična metoda, ki temelji na metodi robnih elementov za izračun tokovnega polja tekočine in Lagrangeovega pristopa za sledenje delcev [13]. Računski algoritem je bil razvit in preizkušen za sledenje delcev pod vplivom hidrodinamičnih in magnetnih sil . Rezultati so pokazali, da je razvit numerični algoritem primeren za reševanja toka razredčene suspenzije in gibanja delcev v magnetno statičnem polju v naprednih naprav za ločevanja delcev. Naprave za ločevanja delcev v obliki ravnega kanala je pokazala slabo učinkovitost izločanja v sredini kanala. To je posledica simetrične porazdelitve magnetnih sil po kanalu, kot tudi simetrični tokovnega polja v kanalu. Prispevek predstavlja nadaljnji razvoj numeričnega algoritma in poglobljeni študiji spremenjene oblike ločevalnega kanala z vstavljenimi periodično ponavljajočimi se pregradami, ki tvorijo v samem kanalu bifurkacije. Cilj spremembe geometrije je ustvariti tokovne motnje in posledično učinkovitejše izločevanje namagnetenih delcev. Taka sprememba geometrije ravnega kanala je pred kratkim pridobila pozornost in je uporabljena v študiji, pretoka krvi v cevi z bifurkacijo [14]. 2. Opis problema Za povečanje učinkovitosti izločevanja delcev (razmerje med številom magnetnih delcev, ki so se prilepili na stene in delci ki so vstopili v izločevalno napravo) iz toka suspenzije smo spremenili geometrijo ravnega kanala kvadratnega preseka [13] v obliko z vstavljenimi trapezoidnimi pregradami (slika 1). Fakulteta za strojništvo, Univerza v Mariboru E-naslov: matej.zadravec@um.si; Tel.: +386-2-220-7783 * Avtor za korespondenco Pomurska obzorja 1/ 2014/ 2 | 87 Matej ZADRAVEC, Matjaž HRIBERŠEK, Jure RAVNIK: LOČEVANJE MAGNETNIH DELCEV IZ… Gostoto magnetnega polja smo spreminjali v območju vrednosti od 0,04T do 0,62T (tabela 2). Suspenzija je bila obravnavana kot razredčena, zato smo zanemarili interakcije med delci. Tabela 2. Podatki koeficient magnetne jakosti Cpm za različne vrednosti Re števila in gostote magnetnega polja B0 . Slika 1. Geometrija ločevalne naprave z vstavljenimi pregradami. Kanal ima kvadratni prečni prerez (0,75x0,75mm) s pregradami dveh različnih višin: večja hbif=h/2=0,0325mm in manjša hbif=h/8=0,0125mm. Dolžina pregrade je enaka višini kanala lbif=H. Razdalja med pregradami je Δlbif=4h=3mm. Magnetno polje je bilo generirano z dvema jeklenima žicama s premerom 0,5mm, ki sta nameščeni na vrhu in na dnu kanala in izpostavljeni zunanjemu homogenemu magnetnemu polju (slika 2). V toku tekočine vode pri temperaturi 20°C so bili razpršene polistirenske kroglice (Micromond Inc, Nemčija) s premerom 1,7μm. Magnetne polistirenske kroglice vsebujejo 12,45% magnetita in 87,55% polistirena [15]. Na vstopu smo predpostavili popolnoma razvit tok in naključno razporeditev delcev po površini vstopa. Zaradi visokega razmerja med višino in dolžino kanala (1:135) bi bili računski časi precej dolgi in smo za reševanje tokovnega polja tekočine območje razdelili na 27 krajših odsekov (slika 1). Segmenti so pozneje pri simulaciji obnašanja delcev združeni s pomočjo periodične robnega pogoja na vstopu in izstopu iz vsakega segmenta (slika 3). Računska mreža obravnavanega problema s pregradami je bila sestavljena iz heksaedrskih elementov in sicer 60x12x12 elementov skupaj torej 186792 računskih vozlišč. Slika 2. Magnetno polje. 3. Numerični model Simulacije so bile izvedene za različne pretoke suspenzije in različne gostote magnetnega polja, ki so bili izbrani na podlagi objavljenih prispevkov [13 in 15] s katerimi smo nato primerjali rezultate. Tok suspenzije je bil laminarni in pretok smo spreminjali od 27ml/h do 162ml/h. Ti pretoki ustrezajo povprečnim hitrostim od u0=1,35cm/s do u0=8cm/s in pripadajočim Reynoldsovim številom od vrednosti 10 do 59,36 (tabela 1). V suspenziji smo upoštevali skupno 10000 polistirenskih delcev. Tabela 1. Podatki za različne hitrosti ter pripadajoče vrednosti brez dimenzijskih števil (Reynolds-ovo in Stokes-ovo število). 88 | Pomurska obzorja 1/ 2014/ 2 Slika 3. Periodični robni pogoji. 4. Rezultati Iz končne porazdelitve izločenih delcev na stenah kanala za primera večje in manjše pregrade v kanalu (slika 4) lahko vidimo, da je učinkovitost izločanja na zgornji in spodnji steni precej večja v primeru manjše pregrade. Matej ZADRAVEC, Matjaž HRIBERŠEK, Jure RAVNIK: LOČEVANJE MAGNETNIH DELCEV IZ… Slika 4. Končni položaj izločenih delcev v kanalu z vstavljenimi pregradami (levo-večja pregrada hbif=h/2 in desno-manjša pregrada hbif=h/8) za primer toka Re=10 in B0=0,62T. Ugotovimo lahko tudi koliko delcev se nam ne izloči na koncu kanala, kar ponovno kaže na slabšo učinkovitost izločanja geometrije z večjo pregrado v primerjavi z manjšo. Tako obnašanje sledi iz dejstva, da je motnja toka v primeru večje pregrade večja, kar vodi do močnejših hidrodinamičnih sil. To še posebej velja v razcepu kanala pri pregradi, kjer bi naj magnetne sile prevladale hidrodinamične sile toka. Učinek je izrazitejši v primeru manjše pregrade saj je le tam hitrost zaradi zožitve manjša. Učinkovitost izločanja delcev ζ na razdalji x (po dolžini kanala) je prikazana na sliki 5 za oba geometrijska primera. Učinkovitost izločanja delcev ζ je definiran kot razmerje med številom magnetnih delcev, ki so se prilepili na stene in delci ki so vstopili v izločevalno napravo. Iz slike 5 izhaja, da je učinkovitost izločanja vzdolž kanala veliko višja v primeru manjše pregrade ¸ Slika 5. Učinkovitost izločanja delcev ζ na razdalji x (po dolžini kanala) za oba geometrijska primera (večja pregrada h bif=h/2 in manjša pregrada hbif=h/8) za primer toka Re=10 in B0=0,62T. Slika 6. Učinkovitost izločanja delcev glede na gostoto magnetnega polja pri hitrosti toka u0=1,35m/s za oba geometrijska primera (večja hbif=h/2 in manjša pregrada hbif=h/8) in kanal brez pregrad [13]. Pomurska obzorja 1/ 2014/ 2 | 89 Matej ZADRAVEC, Matjaž HRIBERŠEK, Jure RAVNIK: LOČEVANJE MAGNETNIH DELCEV IZ… Primerjava učinkovitosti izločanja obeh geometrijskih izvedb naprave s pregradama v kanalu v primerjavi z rezultati izločanja v ravnem kanalu [13] je prikazana na slikah 6 in 7 za različne gostote magnetnega polja in dve različni hitrosti toka u0=1,35m/s (Re=10) in u0=5m/s (Re=37,1). Pri manjšem pretoku Re=10 (slika 6) je tako učinkovitost izločanja največja v primeru geometrijske izvedbe z manjšo pregrado, medtem ko je v primeru višje hitrosti toka Re=37,1 (slika 7) najučinkovitejše izločanje v geometriji brez pregrad [13]. V vseh izračunanih primerih z različnimi hitrostmi toka so najnižje učinkovitost izločanja dosežene z geometrijsko izvedbo, kjer je vključena večja pregrada. Kot je bilo že navedeno, pride do manjše učinkovitosti izločanja pri primeru z večjo pregrado zaradi lokalnega povečanja hitrosti toka nad in pod pregrado. Čeprav so delci izpostavljeni večjim Kelvinovim silam magnetnega polja je v tem območju zadrževalni čas delcev zaradi velikih hitrosti precej manjši, kar sledi že iz podatka o izredno majhnih vrednostih brez dimenzijskega Stokesovega števila delcev (tabela 1) in pomeni, da delci precej dobro sledijo tokovnemu polju. Učinkovitost izločanja je pri vseh simulacijah višja s povečevanjem gostote magnetnega polja in je omejena s točko nasičenja magnetnih lastnosti delcev [13], ki je za nerjavno jekleno žico in magnetne delce B0>0,3T. Slika 7. Učinkovitost izločanja delcev ζ glede na gostoto magnetnega polja pri hitrosti toka u 0=5m/s za oba geometrijska primera (večja hbif=h/2 in manjša pregrada hbif=h/8) in kanal brez pregrad [13]. Slika 8. Učinkovitost izločanja delcev ζ glede na različne hitrosti toka in gostoto magnetnega polja B 0=0,62T za oba geometrijska primera (večja hbif=h/2 in manjša pregrada hbif=h/8) in kanal brez pregrad [13]. 90 | Pomurska obzorja 1/ 2014/ 2 Matej ZADRAVEC, Matjaž HRIBERŠEK, Jure RAVNIK: LOČEVANJE MAGNETNIH DELCEV IZ… Vpliv hitrosti na učinkovitost izločanja je prikazan na sliki 8 za vse geometrijske primere. Povečevanje hitrosti toka zmanjšuje učinkovitosti izločevanja pri vseh geometrijskih izvedbah izločevalne naprave. Iz teh rezultatov je razvidno, da je v primerih manjših hitrosti toka od 2m/s najučinkovitejša izločevalna naprava tista, ki ima vključeno manjšo pregrado. 5. Zaključki S pomočjo 3D simulacij na osnovi metode robnih elementov smo preučevali vpliv geometrijskih lastnosti na učinkovitost izločanja magnetnih delcev v izločevalni napravi (HGMS). Izhajajoč iz prvotne oblike ravnega kanala [13], kjer je učinkovitost izločanja v območjih nizkih magnetnih sil v sredi kanala negativno vplivala na izločanje, je nova zasnova z vstavljenimi ponavljajočimi se trapeznimi pregradami v primeru manjših pregrad povečala učinkovitost izločanja. Rezultati izračunov kažejo, da povečanje velikosti pregrade ne izboljša učinkovitosti izločanja, kot je bila učinkovitost v kanalu brez pregrade. Izvedena raziskava kaže, da ima HGMS naprava s pregrado potencial zlasti v primeru bolj viskoznih tekočin, kjer so vrednosti brez dimenzijskega Reynoldsovega števila Re<10. Reference 1. T. Y. Ying, S. Yiacoumi, C. Tsouris, High-gradient magnetically seeded filtration, Chem Eng Sci 55 (2000) 1101–13. 2. G. D. Moeser, K. A. Roach, W. H. Green, T. A. Hatton, High-gradient magnetic separation of coated magnetic nanoparticles, AIChE J 50 (2004) 2835–48. 3. A. Ditsch, S. Lindenmann, P. E. Laibinis, D. Wang, T. A. Hatton, High-gradient magnetic separation of magnetic nanoclusters, Ind Eng Chem Res 44 (2005) 6824–36. J. Svoboda, V. Ross, Particle capture in the matrix of a magnetic separator, Int J Miner Process 27 (1989) 75–94. 4. 5. G. R., Macroscopic model of particles capture by the elliptic cross-section collector in magnetic separator, J Magn Magn Mater 272–276 (2004) 2348–9. 6. G. Hatch, R. Stelter, Magnetic design considerations for devices and particles used for biological high-gradient magnetic separation (HGMS) systems, Journal of Magnetism and Magnetic Materials 225 (2001) 262–276. 7. H. Tsutsui, C.-M. Ho, Cell separation by non-inertial force fields in microfluidic systems, Mechanics Research Communications 36 (2009) 92–103. 8. T.-Y. Liu, S.-H. Hu, D.-M. Liu, S.-Y. Chen, I.-W. Chen, Biomedical nanoparticle carriers with combined thermal and magnetic responses, Nano Today 4 (1) (2009) 52 – 65. 9. H. Chen, M. Kaminski, X. Liu, Y. Xie, C. Mertz, M. Torno, et al, A novel human detoxification system based on nanoscale bioengineering and magnetic separation techniques., Med Hypotheses 68 (2007) 1071–9. Magnetic Nanoparticles. Synthesis and Rotational Dynamics under Applied Magnetic Fields, Langmuir 21 (2005) 11500–11509. 11. H. Chiriac, A.-E. Moga, G. Iacob, O. C. Mungiu, Amorphous magnetic microspheres for biomedical applications, Journal of Magnetism and Magnetic Materials 293 (1) (2005) 28 – 32. 12. B. D. Plouffe, D. K. Nagesha, R. S. DiPietro, S. Sridhar, D. Heiman, S. K. Murthy, L. H. Lewis, Thermomagnetic determination of fe3o4 magnetic nanoparticle diameters for biomedical applications, Journal of Magnetism and Magnetic Materials 323 (17) (2011) 2310 – 2317. 13. J. Ravnik, M. Hriberˇsek, High gradient magnetic particle separation in viscous flows by 3D BEM, Comput Mech 51 (2013) 465–474. 14. X. Li, A. Popel, G. Karniadakis, Bloodplasma separation in y-shaped bifurcating microfluidic channels, Physical Biology 9 (2012) 1–12. 15. H. Chen, M. D. Kaminski, A. J. Rosengart, 2D modeling and preliminary in vitro investigation of a prototype high gradient magnetic separator for biomedical applications, Medical Engineering and Physics 30 (2008) 1–8. 10. H. Singh, P. Laibinis, T. Hatton, Rigid, Superparamagnetic Chains of Permanently Linked Beads Coated with Pomurska obzorja 1/ 2014/ 2 | 91 Biotehnika Milan Šernek1,*, Maksimilijan Mravljak1 Adhezija pri lepljenju lesa z jeklom POVZETEK Pri montaži lepljenih lameliranih nosilcev v želeno konstrukcijo se uporabljajo tudi nove vrste vezi - jeklene palice z metričnim navojem, ki jih s pomočjo lepila usidrajo oziroma vlepijo v leseni nosilec. Imenujemo jih »glued-in rods« ali vezi iz vlepljenih jeklenih palic (VIVJP). Za izvedbo kakovostnih VIVJP je ključnega pomena razumevanje kohezije in adhezije pri lepljenju lesa in jekla, za kar je potrebno poznati reološke lastnosti lepila ter relevantne lastnosti obeh lepljencev (les in jeklo). Glavni dejavniki, ki vplivajo na trdnost vezi, so geometrija (globina in premer) in morfologija izvrtine, vlažnost lesa, anatomska smer lesa in potek lesnih vlaken ter debelina in površina lepilnega spoja. Ključne besede: adhezija, epoksidno lepilo, jeklo, lepljeni lamelirani nosilci, les. 1. Uvod Lepljeni lamelirani nosilci »glued laminated timber« (glulam) so namenjeni gradnji objektov z večjimi razponi, kot so športne dvorane, pokriti bazeni, paviljoni, garažne hiše, industrijske hale in podobno. Dimenzije nosilcev so omejene s proizvodnimi zmogljivostmi, določene omejitve pa obstajajo tudi glede transporta od tovarne do gradbišča. Pogosto je potrebno lepljene lamelirane konstrukcijske elemente spojiti med sabo, za kar se večinoma uporabljajo vezi na osnovi vijakov z maticami ali mozniki (Mravljak in Šernek, 2009). Pri večjih konstrukcijah je lahko število vijakov ali moznikov za posamezen spoj zelo veliko, kar predstavlja dodatno težavo pri sestavi. Takšne vezi z estetskega vidika niso najbolj zaželene (Kuhlmann in sod., 2001), zato se za spajanje lesa in lesnih kompozitov v kompleksnejšo konstrukcijo uporabljajo nove vrste vezi, ki so se zadnjem času uveljavile v praksi. Take konstrukcijske vezi, ki se uporabljajo predvsem pri montaži lepljenih nosilcev, predstavljajo jeklene palice z metričnim navojem, ki jih s pomočjo lepila usidrajo oziroma lepijo v les. Imenujemo jih »glued-in rods« ali vezi iz vlepljenih jeklenih palic (VIVJP). Za izvedbo kakovostnih VIVJP je ključnega pomena razumevanje kohezije in adhezije pri lepljenju lesa in jekla, za kar je potrebno poznati reološke lastnosti lepila (Mravljak in Šernek, 2011) ter relevantne lastnosti lepljencev (les in jeklo). Pomemben vpliv na obnašanje teh vezi imajo tudi geometrija 1Univerza v Ljubljani, Biotehniška fakulteta, Oddelek za lesarstvo, Rožna dolina C. VIII/34, 1000 Ljubljana E-naslov: milan.sernek@bf.uni-lj.si * Avtor za korespondenco; Tel.: +1-320-3623; Fax: +1-257-2297 92 | Pomurska obzorja 1/ 2014/ 2 (globina in premer) in morfologija izvrtine, ki je odvisna od načina izdelave, anatomska smer lesa in potek lesnih vlaken, vlažnost lesa ter drugi dejavniki, kot sta debelina lepilnega sloja oz. spoja in način obremenitve. Cilj predstavljene raziskave je bil proučiti vpliv omenjenih dejavnikov na trdnost vezi in s tem zagotoviti varno uporabo VIVJP. 2. Materiali in metode 2.1. Izbira materialov V raziskavi smo uporabili tri komercialna epoksidna lepila, primerna za lepljenje lesa in jekla: 1. 2. 3. ERGO 7211 EPOXY 2K 50 ml (7211). EPOX 210 A+B (EPOX 210). NEOSTIK EP 101 (EP 101). Za izdelavo testnih preskušancev smo uporabili smrekov les (Picea abies L.), jeklene ploščice in jeklene palice z metričnim navojem. 2.2. Proučevanje reoloških lastnosti lepil V začetnem delu raziskave smo raziskali reološke lastnosti treh epoksidnih lepil med utrjevanjem pri različnih temperaturah in nato izbrali najustreznejše lepilo za nadaljnje raziskave lepljenja lesa in jekla. Izvedli smo oscilatorni test z reometom ARES G2 (slika 1). Milan ŠERNEK, Maksimilijan MRAVLJAK: ADHEZIJA PRI LEPLJENJU LESA Z JEKLOM Slika 1. Zasnova oscilatornega testa z reometrom (levo) in reometer ARES G2 (desno). 2.3. Preskušanje strižne trdnosti lepilnega spoja Za proučevanje trdnosti vezi med lesom in jeklom smo izdelali strižne preskušance s preklopom, in sicer tako, da je bil en del preskušanca iz lesa, drugi del pa iz jekla (slika 2). Oba dela preskušanca smo med seboj zlepili z epoksidnim lepilom ERGO 7211. Leseni del preskušanca (lamela dimenzij 150 mm x 20 mm x 5 mm) je bil izdelan iz smrekovine z ravnovesno vlažnostjo 12 %. Izdelali smo tri vrste lamel glede na orientacijo lesnih vlaken: 1. Lamele z radialno strukturo 2. Lamele s tangencialno strukturo 3. Lamele z radialno/tangencialno strukturo Jekleni del preskušanca je bil izdelan iz jekla z nominalno natezno trdnostjo 1000 N/mm2. Strižni test smo opravili v skladu s standardnim postopkom z univerzalnim testirnim strojem Zwick Z100. 2.4. Preskušanje VIVJP V zaključnem delu raziskave smo izdelali VIVJP ter ugotavljali vpliv debeline lepilnega sloja in sidriščne dolžine na trdnost vezi. Leseni del VIVJP smo izdelali iz smrekovega lesa in pripravili šest ponovitev za vsako od naslednjih kombinacij preskušancev: 1. Tri različne debeline lepilnega sloja okoli jeklene palice: 0 mm, 0,5 mm in 1 mm. 2. Dva različna premera palic: 10 mm in 12 mm. 3. Dve sidriščni dolžini: 100 mm in 150 mm. Presek preskušancev smo določili glede na premer vlepljene palice pri čemer smo upoštevali varnostno razdaljo od roba preskušanca, ki je potrebna, da ne pride do nezaželenega prečnega pokanja po lesu. Jekleni del preskušanca smo izdelali iz palic z metričnim navojem premera 10 mm in 12 mm. Za preskušance z globino izvrtine 100 mm smo pripravili palice dolžine 150 mm, za preskušance z globino izvrtine 150 mm pa palice dolžine 200 mm. Tako smo zagotovili 50 mm jeklene palice izven lesenega dela, potrebne za učinkovito vpetje v čeljust trgalnega stroja (slika 3). Palice smo temeljito očistili in razmastili ter tako preprečili vnos nečistoč v lepilo oz. spoj. Lepili smo z epoksidnim lepilom ERGO 7211, ki smo ga nanesli v izvrtino, vanjo vstavili jekleno palico in pustili, da je lepilo utrdilo. Slika 2. Sestava strižnih preskušancev s preklopom s pomočjo šablone (levo) in zlepljen preskušanec po porušitvi (desno). Pomurska obzorja 1/ 2014/ 2 | 93 Milan ŠERNEK, Maksimilijan MRAVLJAK: ADHEZIJA PRI LEPLJENJU LESA Z JEKLOM Slika 3. Vpetje VIVJP v univerzalni testirni stroj Zwick Z100 (levo), detajl vpetja vijaka z metričnim navojem (v sredini) in porušitev VIVJP preskušanca po lesnih vlaknih (desno). 3. Rezultati in razprava Preglednica 1. Časi želiranja (tžel) in zamreženja (tzam) proučevanih lepil pri različnih temperaturah (T). 3.1. Utrjevanje lepil Lepilo S pomočjo reometrije smo spremljali proces utrjevanja lepil in določili njihove karakteristične lastnosti. Profil utrjevanja proučevanih epoksidnih lepil pri temperaturi 30 °C je prikazan na sliki 4. Na osnovi sečišča krivulj G' in G'' (točka želiranja) smo ugotovili, da je najhitreje utrdilo lepilo ERGO 7211, nato EPOX 210 in zatem EP 101. EP 101 EPOX 210 ERGO 7211 T (°C) tžel (min) tzam (min) tžel (min) tzam (min) tžel (min) tzam (min) 30 499 565 345 527 238 297 40 186 365 151 321 139 174 50 91 277 69 272 64 103 60 44 334 38 457 30 57 80 8 103 11 1183 6 33 3.2. Strižna trdnost lepilnega spoja 1,E+07 G'-EP 101 G'-7211 9,E+06 G'-EPOX 210 Ugotovili smo, da orientacija lesnih vlaken v lesenem delu preskušanca ni bistveno vplivala na strižno trdnost lepilnega spoja (slika 5). Delež loma po lesnih vlaknih je bil v vseh primerih zelo visok (slika 6), kar pomeni, da je bil lepilni spoj »močnejši«, in da ugotovljena trdnost pravzaprav predstavlja strižno trdnost smrekovine. S praktičnega vidika to pomeni, da je bilo izbrano epoksidno lepilo kakovostno in primerno za lepljenje lesa in jekla. 8,E+06 G', G'' (Pa) 7,E+06 6,E+06 5,E+06 4,E+06 3,E+06 G''-EP 101 G''-7211 2,E+06 G''-EPOX 210 1,E+06 0,E+00 0 120 240 360 480 600 720 840 960 1080 1200 1320 1440 1560 Čas (min) Slika 4. Prikaz profilov utrjevanja proučevanih epoksidnih lepil pri 30 °C. Temperatura utrjevanja je imela bistveni vpliv na čas želiranja in zamreženja lepil (preglednica 1). Z araščajočo temperaturo se je povečala hitrost utrjevanja lepil in s tem skrajšal čas za želiranje in zamreženje. Slika 5. Strižna trdnost lepilnega spoja glede na orientacijo lesenega dela preskušanca. 94 | Pomurska obzorja 1/ 2014/ 2 Milan ŠERNEK, Maksimilijan MRAVLJAK: ADHEZIJA PRI LEPLJENJU LESA Z JEKLOM RILEM - Symposium on Joints in Timber Structures, 2001, Stuttgart, Germany: 343-352. 2. Mravljak M., Šernek M. Vpliv osnovnih dejavnikov na kratkotrajno obnašanje vezi iz vlepljenih jeklenih palic. Les, 2009, 61, 7-8. 3. Mravljak M., Šernek M. The influence of curing temperature on the rheological properties of epoxy adhesives. Drvna industrija, 2011, 62(1), 19-25. Slika 6. Visok delež loma po lesu preskušancev s preklopom 3.3. Trdnost VIVJP Maksimalna obremenitev potrebna za porušitev VIVJP je naraščala z sidriščno dolžino jeklene palice in njenim premerom (preglednica 2). To je bilo pričakovano, saj se je pri tem povečala lepilna površina. Na trdnost vezi je vplivala tudi debelina lepilnega sloja in sicer so bili v večini primerov najboljši rezultati doseženi pri srednji debelini (0,5 mm). Preglednica 2. Maksimalna obremenitev (kN) VIVJP glede na debelino lepilnega sloja, premer palice in sidriščno dolžino. Debelina lepilnega Premer palice  in sidriščna dolžina sloja (mm) (mm) Ø 10 Ø 12 100 150 100 150 0,0 21,9 23,5 23,4 29,7 0,5 22,9 23,1 29,6 33,3 1,0 20,2 18,8 25,6 29,7 4. Sklepi Ugotovili smo, da lahko z merjenjem reoloških lastnosti zelo podrobno spremljamo proces utrjevanja epoksidnih lepil ter določimo čas želiranja in zamreženja. Z naraščajočo temperaturo utrjevanja sta se omenjena časa bistveno skrajšala. Najhitreje je utrdilo lepilo ERGO 7211, nato EPOX 210 in najkasneje EP 101. Orientacija lesnih vlaken ni vplivala na strižno trdnost preskušancev s preklopom. Lom preskušancev je potekal po lesnih vlaknih. Največjo obremenitev so prenesle VIVJP z daljšo (150 mm) sidriščno dolžino, večjim (12 mm) premerom palice in s srednjo (0,5 mm) debelino lepilnega sloja. Zahvala Raziskava je potekala v okviru raziskovalnega programa P40015 »Les in lignocelulozni kompoziti«, ki ga financira Javna agencija za raziskovalno dejavnost RS. Literatura 1. Kuhlmann U., Aicher S., Lippert P. Rigid frame corners with glued-in rods. V: Proceedings of the International Pomurska obzorja 1/ 2014/ 2 | 95