["itemContainer",{"xmlns:xsi":"http://www.w3.org/2001/XMLSchema-instance","xsi:schemaLocation":"http://omeka.org/schemas/omeka-xml/v5 http://omeka.org/schemas/omeka-xml/v5/omeka-xml-5-0.xsd","uri":"http://romeka.rgf.rs/items/browse?tags=deterministic+chaos&output=omeka-json","accessDate":"2020-01-21T09:29:20+01:00"},["miscellaneousContainer",["pagination",["pageNumber","1"],["perPage","500"],["totalResults","1"]]],["item",{"itemId":"762","public":"1","featured":"0"},["fileContainer",["file",{"fileId":"917"},["src","http://romeka.rgf.rs/files/original/Doktorske_disertacije/DD_Kostic_Srdjan/DD_Kostic_Srdjan.1.pdf"],["authentication","46d2d4ea7aa91279d346fa00fc4de4f1"],["elementSetContainer",["elementSet",{"elementSetId":"5"},["name","PDF Text"],["description"],["elementContainer",["element",{"elementId":"133"},["name","Text"],["description"],["elementTextContainer",["elementText",{"elementTextId":"17035"},["text","UNIVERZITET U BEOGRADU RUDARSKO-GEOLOŠKI FAKULTET\r\nSrđan D. Kostić\r\nNELINEARNO DINAMICKO MODELOVANJE POTRESA IZAZVANIH PROMENOM NAPONSKOG STANJA PRI IZRADI HORIZONTALNIH PODZEMNIH PROSTORIJA\r\nDoktorska disertacija\r\nBeograd, 2013.g.\r\nMentor:\r\nDr Nebojša Vasović, vanr. prof.,\r\nUža naučna oblast: Mehanika\r\nUniverzitet u Beogradu Rudarsko-geološki fakultet\r\nČlanovi komisije:\r\n1. Dr Slobodan Trajković, red. prof.\r\nUža naučna oblast: Eksploatacija čvrstih mineralnih sirovina i mehanika stena Univerzitet u Beogradu Rudarsko-geološki fakultet\r\n2. Dr Nebojša Gojković, red. prof.\r\nNaučna oblast: Eksploatacija čvrstih mineralnih sirovina i mehanika stena Univerzitet u Beogradu Rudarsko-geološki fakultet\r\n3. Dr Dragutin Jevremović, red. prof.\r\nNaučna oblast: Mehanika tla i fundiranje, mehanika stena i inženjerska geologija Univerzitet u Beogradu Rudarsko-geološki fakultet\r\n4. Dr Nikola Burić, naučni savetnik\r\nUža naučna oblast: Kvantna i matematička fizika Univerzitet u Beogradu Institut za fiziku\r\nDatum odbrane:\r\nNELINEARNO DINAMIČKO MODELOVANJE POTRESA IZAZVANIH PROMENOM NAPONSKOG STANJA PRI IZRADI HORIZONTALNIH PODZEMNIH PROSTORIJA\r\nREZIME\r\nU ovoj disertaciji, cilj istraživanja bio je određivanje dinamičkih svojstava potresa izazvanih promenom naponskog stanja pri izradi horizontalnih podzemnih prostorija. Analiza potresa izvedena je sa dva aspekta - u prvom delu istraživanja izučavane su vremenske serije oscilovanja tla za vreme potresa, a u drugom delu vršena je analiza stabilnosti fenomenološkog Baridž-Knopof modela bloka sa oprugom. Vremenske serije oscilovanja tla proučavane su primenom postupka analize surogat podataka, kao i korišćenjem metoda analize nelinearnih vremenskih serija (primena teoreme o razvijanju skalarne vremenske serije). Analiza vremenskih serija izvedena je za potrese izazvane reaktiviranjem kretanja duž raseda usled napredovanja rudarskih radova, registrovanih na tri lokacije (rudnika) u Poljskoj: rudnici bakra Rudna i Legnic-Glogov, i rudnik uglja Gornja Silezija. Rezultati izvedene analize pokazali su da registrovane vremenske serije oscilovanja tla u pravcu sever-jug i u vertikalnom pravcu pripadaju nelinearnim vremenskim serijama, dok oscilovanja tla u pravcu istok-zapad pripadaju klasi stacionarnih stohastičkih procesa sa Gausovom raspodelom stohastičkog dela, koji mogu biti modifikovani nepoznatom nelinearnom funkcijom. Međutim, iako su oscilovanja tla u pravcu sever-jug i u vertikalnom pravcu nelinearna, primenom teoreme o razvijanju potvrđena je njihova stohastička priroda, i to malom vrednošću determinističkog faktora K (<1) i niskom greškom unakrsnog predviđanja kod testa stacionarnosti. U drugom delu istraživanja, vršena je analiza stabilnosti kretanja nekoliko varijanti Baridž-Knopof modela, i to: Madariaga modela jednog bloka sa Diterih-Ruina zakonom trenja zavisnim od brzine i stanja, potom Karlson-Langer modela sa jednim i dva bloka i sa zakonom trenja zavisnim samo od brzine, kao i Beker modela sa jednim blokom, Diterih-Ruina zakonom trenja zavisnim od brzine i stanja i sa dve promenljive stanja. Analiza stabilnosti izvedena je numeričkim i analitičkim putem. Numerički je izvršena analiza stabilnosti Karlson-Langer i Beker modela, posmatranjem ponašanja rešenja u blizini ravnotežnog stanja, dok je u slučaja Madariaga modela, bifurkaciona analiza izvedena analitički (standardna lokalna bifurkaciona analiza), a rezultati su potvrđeni numerički u programskom paketu DDE-BIFTOOL. Analiza dinamike ovih modela vršena je u uslovima promene uvedenih parametara i to: parametra jačine trenja C u Karlson-Langer modelu, kao i parametra vremenskog kašnjenja T u izrazu za trenje u Karslon-Langer i Madariaga modelu. Pored toga, dinamika modela je osmatrana i za periodičnu perturbaciju pojedinih parametara sistema, u Madariaga i Bekerovom modelu. Opravdanje za uvođenje novih parametara, kao i za njihovu oscilatornu promenu može se pronaći u sledećem. Parametar vremenskog kašnjenja T na kvalitativan način modeluje ,,zakasneli'' uticaj napredovanja iskopa podzemne prostorije i promene naponskog stanja na reaktivaciju kretanja duž raseda. Takođe, ovim parametrom se opisuje uticaj stanja površi kontakta između bloka i hrapave podloge na dinamiku kretanja. Parametar jačine trenja C vezan je za inherentna svojstva rasedne površi - debljinu rasedne zone, njenu dubinu, mineralni sastav i svojstva minerala, i sl. Periodična perturbacija parametara odražava uticaj napredovanja rudarskih radova na reaktiviranje kretanja duž raseda i pojavu seizmičke aktivnosti (vibracija rudarskih mašina, efekti miniranja, i sl.). Izvedena analiza je pokazala da, sa uvedenim vremenskim kašnjenjem, izučavani modeli pokazuju determinističko haotično ponašanje, koje predstavlja tipičnu (generičku) pojavu u Madariaga modelu, dok se, u Karlson-Langer modelu javlja samo kao tranzijentna pojava. Neophodno je naglasiti da uvođenje vremenskog kašnjenja, osim indukovanja kompleksnijeg ponašanja, može uzrokovati i stabilizaciju kretanja bloka, i njegov povratak u ravnotežno stanje (oscilatorna smrt). U pogledu ponašanja sistema pri periodičnoj perturbaciji parametara, i u Madariaga i u Beker modelu, ovakva promena izaziva pojavu determinističkog haotičnog ponašanja, i to za mnogo manje vrednosti parametara u odnosu na vrednosti dobijene prethodnim istraživanjima. Pri tome, generisanje kompleksnog ponašanja uočeno je pri oscilacijama malih amplituda, što odgovara mogućem uticaju napredovanja rudarskih radova (vibracija rudarskih mašina, i sl.). Takođe, kompleksna dinamika sistema osmatrana je za vrednosti frekvencija oscilacija koje su bliske frekvenciji oscilacija bloka u neperturboavnom stanju. Drugim rečima, do reaktiviranja kretanja duž raseda može doći i usled oscilacija čije su amplitude male, a frekvencije približno iste sa frekvencijom aseizmičkog kretanja duž raseda. Opšte posmatrano, izvedeno istraživanje je pokazalo da je vremenski interval između iskopa podzemne prostorije, promene naponskog stanja i reaktiviranja kretanja duž raseda od ključne važnosti za karakter dinamike kretanja duž raseda. Neophodno je naglasiti da se prividni nesklad između stohastičke prirode oscilacija tla za vreme potresa i determinističkog haotičnog modela mehanizma nastanka potresa može objasniti činjenicom da procesi nastanka potresa i kretanja talasa kroz Zemljinu koru predstavljaju disipativne procese. Dakle, pri nastajanju potresa, jedan deo energije se troši na savlađivanje trenja duž raseda, frikciono zagrevanje, i dr., dok samo mali deo biva emitovan u vidu seizmičkih talasa. Potom se energija tih talasa, pre registrovanja na površini terena, rasipa na svom putu kroz Zemjinu koru, usled složenih tektonskih odnosa, geološke građe, morfološki izraženih struktura, distorzije talasa na tom putu, i sl. Izvedeno istraživanje predstavlja čvrstu osnovu za nastavak daljih istraživanja, kroz primenu savremenih tehnologija i povećanjem broja osmatračkih stanica, sa jedne strane, kao i kroz uspešniju analizu dobijenih podataka i ponašanja predloženih fenomenoloških modela, sa druge strane.\r\nKLJUČNE REČI: potres, naponsko stanje, podzemne prostorije, bifurkacije, vremenska serija, deterministički haos, stohastičnost, perturbacije, vremensko kašnjenje, jačina trenja\r\nNAUČNA OBLAST: Rudarsko inženjerstvo\r\nUŽA NAUČNA OBLAST: Eksploatacija čvrstih mineralnih sirovina i mehanika stena\r\nUDK broj: 512.643:519.245/.856/.87:534:64:550.8:616.232\r\n622.1//.261.2/.268.8/.273/.33/.83\r\n624.075/.121.5/.131.55\r\n(043.3)\r\n1. UVOD\r\nRudarsko inženjerstvo, kao naučna oblast, predstavlja sintezu različitih naučnih disciplina, pre svega prirodnih i tehničkih, tako da se tokom poslednjih godina sve više nameće potreba za primenom savremenih matematičkih metoda i tehnika, koje mogu doprineti boljem razumevanju geoloških pojava i procesa, među kojima i rudarski izazvanih potresa, kao najvažnijih prirodnih faktora koji utiču na uslove izgradnje i eksploatacije rudarskih objekata. Tradicionalni pristupi, koji podrazumevaju primenu empirijskih formula, ne daju zadovoljavajuće rezultate, imajući u vidu činjenicu da su pomenute formule izvedene za pojedine vrste stenskih masa, i, samim tim, imaju ograničenu primenu u drugim geološkim sredinama. U tom smislu, prateći, pre svega svetske trendove u oblasti izučavanja zemljotresa, došlo se do saznanja da bi tehnike i metode izučavanja ovih pojava mogle biti primenjene i na potrese izazvane promenom naponskog stanja pri izradi podzemnih prostorija, što bi dovelo do boljeg razumevanja samog mehanizma nastanka procesa, kao i dinamičkih karakteristika njihovih pojava, i na taj način dalo osnovu za racionalnije projektovanje uslova i tehnologije iskopa, načina podgrađivanja, vrste i karakteristika podgrade, kao i adekvatnih preventivnih i sanacionih mera. Do sada obavljena ispitivanja u svetu, na primeru definisanja mehanizma prirodnih zemljotresa, kao i u drugim oblastima geološke i rudarske nauke, samo potvrđuju hipotezu da je neophodno primeniti moderne matematičke metode i u oblasti rudarski generisanih potresa, među kojima i metode nelinearne dinamike.\r\nZa potrebe eksploatacije mineralnih sirovina, ili tokom izvođenja podzemnih radova za potrebe izrade podzemnih skladišta, skloništa, i drugih vrsta podzemnih objekata, naročito na većim dubinama, remeti se postojeće naponsko stanje u okolnoj stenskoj masi. U slučaju postojanja raseda ili sličnog mehaničkog (ili litološkog) diskontinuiteta u širem području izvođenja rudarskih radova, realno je očekivati da promena napona uzrokuje kretanje duž postojećih diskontinuiteta, čime bi se stvorili uslovi za nastanak potresa. Intenzitet ovih potresa postaje veći i njihovi efekti značajniji sa daljim napredovanjem iskopa i povećanjem dubine, čime se povećavaju i njihove posledice kako po sam proces rada, tako i na materijalna dobra na površini terena. Ova veza između napredovanja rudarske aktivnosti i pojave potresa utvrđena je u mnogim rudarskim oblastima u svetu. Tako, na primer, očigledna je veza između količine iskopanog uglja i registrovane energije seizmičke aktivnosti u Rurskoj rudarskoj oblasti u Nemačkoj (slika 1-1).\r\nSlika 1-1. Seizmička energija i zapremina otkopanog uglja po mesecima za Rursku oblast u Nemačkoj [Fritschen, 2010; modifikovano].\r\nOvakva uzročno-posledična veza između tehnogenih potresa i napredovanja rudarskih radova od posebnog je značaja za teritoriju naše zemlje, s obzirom na to da se pojedina ležišta mineralnih sirovina, u kojima se vrši ili se predviđa podzemna eksploatacija, nalaze u seizmički aktivnim zonama (Baljevac, Jarando, Bela stena, Rudnik, Suvo rudište, Zajača, Krupanj, i dr.). U ovim zonama realno je očekivati da usled povećanog obima eksploatacije mineralnih sirovina i ubrzanog napredovanja radova dođe do pojačane frekvencije potresa, ili do povećanja njihovih magnituda, a što upravo može biti uslovljeno promenom naponskog stanja u stenskoj masi i reaktiviranjem kretanja duž postojećih ruptura u Zemljinoj kori. Neophodno je naglasiti da do sada u našoj zemlji nisu zabeleženi potresi ovog tipa. Jedino su u rudniku uglja u Zenici, u jami ,,Raspotočje'', u dva navrata, rokom 2010.g. i 2012.g. zabeleženi gorski udari, izazvani naglim oslobađanjem velike količine stenske mase iz čela iskopa. S druge strane, postoji veliki broj primera ovakvih potresa u svetu, naročito u zemljama sa razvijenom rudarskom aktivnošću, poput Poljske, Nemačke, SAD-a i Južne Afrike, gde usled intenzivne rudarske aktivnosti često dolazi do pojave potresa različitih magnituda. U Prilogu br.1 dat je pregled registrovanih rudarski generisanih potresa, po kontinentima, državama, rudarskoj oblasti, vrsti sirovine koja se eksploatiše, mehanizmu nastanka i rasponu magnituda.\r\nIzučavanje rudarski izazvanih potresa takođe je značajno i zbog velikih materijalnih šteta koje uzrokuje ovaj tip potresa. Na slici 1-2 dat je primer šteta, izraženih u brojevima dana bez rudarske aktivnosti u rudniku zlata Klerksdorp u Južnoj Africi, u zavisnosti od magnitude i vrste potresa. Kao što se može videti, značajan udeo rudarski izazvanih potresa predstavljaju potresi nastali reaktiviranjem kretanja duž raseda.\r\nSlika 1-2. Veza između pričinjene štete, izražene u brojevima dana bez rudarske aktivnosti, i magnitude različitih tipova potresa u rudniku zlata Klerksdorp u Južnoj Africi [Gay, 1993; modifikovano].\r\nU ovoj disertaciji, predmet istraživanja predstavljaju potresi izazvani promenom naponskog stanja pri izradi horizontalnih podzemnih prostorija. Naglasak na promeni naponskog stanja kao uzročnika pojave potresa proizilazi iz činjenice da i deformacija stenske mase dovodi do pojave potresa. Međutim, potresi ovakvog tipa su manje magnitude i mogu se izbeći pravilnim izborom tehnologije iskopa, definisanjem odgovarajućih parametara fizičko-mehaničkih svojstava stenske mase, i dr. Naime, ukoliko usvojimo podelu rudarski izazvanih potresa predloženu u radu [Hasegawa, i dr., 1989], izdvajaju se šest osnovnih modela nastanka potresa indukovanih podzemnom eksploatacijom: potresi izazvani provaljivanjem krovine podzemnih prostorija, lomom zaštitnih stubova, prekoračenjem čvrstoće na zatezanje povlatne stenske mase, kao i reaktiviranjem kretanja duž normalnog, reversnog i transkurentnog raseda. Ovi modeli se obično razmatraju pri opštoj analizi rudarski izazvane seizmičnosti [Gibowicz, Kijiko, 1994]. Pri tome, pojava prva tri tipa potresa može biti izbegnuta, ili njihova veličina može biti u značajnoj meri redukovana, primenom odgovarajućeg načina eksploatacije, pravilnim odabirom vrste podgrade, pre svega u zavisnosti od vrste i fizičko-mehaničkih svojstava stenske mase. S druge strane, potresi izazvani reaktiviranjem kretanja duž raseda ne mogu se kontrolisati, što ih čini jednim od najznačajnijih geoloških hazarda u područjima sa podzemnom eksploatacijom mineralnih sirovina. Ovakva kompleksnost rudarski generisanih potresa duž postojećih raseda proizilazi i iz njihove ,,haotične'' dinamike mehanizma nastanka potresa, kao i iregularnog kretanja (oscilovanja) tla za vreme potresa.\r\nNeophodno je naglasiti da je kod svih predloženih modela nastanka rudarski izazvane seizmičnosti u radu [Hasegawa, i dr., 1989] razmatrana promena naponskog stanja izazvana izradom horizontalnih podzemnih prostorija za potrebe daljeg napredovanja rudarskih radova, zbog čega je ovaj vid rudarskih radova i predmet doktorske disertacije. U daljem tekstu, ovi potresi biće nazivani ,,rudarski izazvanim potresima'', pod kojima će se podrazumevati potresi nastali promenom naponskog stanja pri izradi horizontalnih podzemnih prostorija.\r\nIzučavanje potresa, u ovoj disertaciji, izvodi se na dva načina: analizom dinamike registrovanih vremenskih serija potresa, i njihovog mehanizma nastanka, oko nepodgrađene podzemne prostorije kružnog poprečnog preseka i u neposrednoj blizini čela otkopa. Cilj istraživanja je dvostruk: (1) određivanje područja u parametarskom prostoru sistema sa različitim dinamičkim ponašanjem, uz pretpostavku periodičnih promena posmatranih parametara i još vremenskog kašnjenja između promene naponskog stanja i reaktiviranja kretanja duž raseda; (2) ispitivanje dinamike oscilovanja tla za vreme potresa.\r\nMetode koje su korišćene u istraživanju već su potvrđene u oblasti izučavanja zemljotresa, s obzirom na to da do nastanka ovih potresa dolazi kretanjem duž postojećeg raseda, a što predstavlja mehanizam koji odgovara tektonskim zemljotresima [Gibowicz, Kijiko, 1994]. Ispitivanje je izvedeno korišćenjem metoda nelinearne dinamike, i to u dve faze. U prvoj fazi, izvršena je analiza nelinearnih vremenskih serija registrovanih potresa na nekoliko izabranih lokacija (rudnika). Analiza je izvedena primenom modernih metoda analize nelinearnih vremenskih serija, kao i postupka ispitivanja surogat podataka. Potom, u drugoj fazi, za potrebe modelovanja mehanizma nastanka rasednih potresa izazvanih napredovanjem rudarskih radova u Zemljinoj kori, pristupilo se analitičkom i numeričkom ispitivanju dinamike modela nastanka potresa. Međutim, neophodno je naglasiti da još uvek ne postoji konsensus oko jedinstvenog modela koji opisuje kretanje duž seizmogenog raseda, a samim tim i nastanak potresa. Stoga, vrlo je korisno izučavati različite modele, razmatrajući pri tom uticaje različitih parametara sistema na dinamiku ispitivanih modela. U ovom istraživanju polazi se od Baridž-Knopof modela bloka sa oprugom1 [Burridge, Knopoff, 1967], čija se dinamika definiše pomoću tri različita sistema diferencijalnih jednačina - prvi predložen u radu [Erickson, i dr. 2008], drugi predložen u radu [Carlson, Langer, 1989], odnosno u radu [Vieira, 1999] i treći predložen u radu [Becker, 2000]. Cilj izvođenja analize jeste utvrđivanje karaktera procesa koji dovodi do pojave potresa. Drugim rečima, neophodno je utvrditi da li je proces koji izaziva pojavu potresa determinističke ili stohastičke prirode. Ukoliko se potvrdi deterministička priroda procesa, onda je moguće nastanak potresa modelovati determinističkim jednačinama, koje su vrlo primenljive u inženjerskoj praksi. S druge strane, u slučaju da analiza dinamičkog ponašanja ukaže na prisustvo determinističkog haosa, to bi direktno ukazivalo na postojanje mogućnosti\r\nkratkoročne predvidivosti procesa za određene vrednosti parametara sistema, sa kvalitativnim korelacijama između teorijsko ispitivanih vrednosti i realno osmatranih podataka.\r\nAnaliza dinamike kretanja modela nastanka potresa izvedena je za tri varijante BK modela, i to:\r\n1) sa jednim i dva bloka, i zakonom trenja samo u funkciji brzine bloka (Karlson- Langer sistem);\r\n2) sa jednim blokom i Diterih-Ruina zakonom trenja [Dieterich, 1979; Ruina, 1983] u funkciji brzine bloka i stanja kontaktne površi (Madariaga sistem);\r\n3) sa jednim blokom i Diterih-Ruina zakonom trenja [Becker, 2000], i dve promenljive stanja (Bekerov sistem).\r\nIzučavanju ovih modela pristupilo se na sledeći na čin:\r\na) uvodeći dva nova parametra:\r\n- vremensko kašnjenje T u izraz za trenje u prvom i drugom modelu;\r\n- parametar jačine trenja c u izraz za trenje u prvom modelu.\r\nb) periodičnom perturbacijom parametara u prvom i trećem modelu. Razlozi za uvođenje vremenskog kašnjenja T su višestruki:\r\n1) vremensko kašnjenje T kvalitativno odražava zakasneli uticaj napredovanja rudarskih radova na reaktiviranje kretanja duž raseda. Naime, promenom naponskog stanja napredovanjem čela iskopa (i time uklanjanjem dela stenske mase), dolazi do promene naponskog stanja i duž rasedne zone. Da bi došlo do kretanja duž raseda, a time i pojave nestabilnosti, potrebno je da, u određenom vremenskom intervalu, promena napona postane veća od koeficijenta trenja duž raseda. Upravo ovaj vremenski interval, ,,hijatus'', između vremena iskopa i vremena potresa kvalitativno se može modelovati uvođenjem vremenskog kašnjenja T;\r\n2) uvođenjem T u izraz za trenje kod modela sa jednim blokom, javljaju se tipovi dinamičkog ponašanja osmatrani u modelima sa više blokova. Na taj način, broj blokova modela se redukuje, bez kvalitativnih promena dinamičkog ponašanja;\r\n3) povećanje statičkog koeficijenta trenja sa vremenom (a time i čvrstoće na smicanje), od ključne je važnosti za pojavu nestabilnosti trenja, a, time i za pojavu potresa.\r\nOpšte posmatrano, čvrstoća na smicanje (koeficijent trenja) predstavlja dinamičku promenljivu, koja je suštinski povezana sa brzim procesima tokom kretanja blokova, i sporim procesima tokom stanja mirovanja blokova [Ben- David, i dr., 2010]. Ovaj spori proces ,,zaceljenja raseda'', koji započinje po prestanku kretanja bloka, obično se naziva ,,starenjem'' raseda, i odvija se tokom vremena koje je određeno ,,lokalnom memorijom'' sistema o efektivnom kontaktu između bloka i hrapave površi. Drugim rečima, proces ,,starenja'' raseda, koji podrazumeva obnavljanje kontakta između bloka i hrapave podloge, kao i logaritamski porast koeficijenta statičkog trenja u dužem vremenskom intervalu, primarno je određen dužinom i karakterom kontakta između bloka u kretanju i hrapave površi donje ploče u prethodnoj brzoj fazi kretanja. Vremensko kašnjenje T, kao dopuna članu za trenje, opisuje ovaj složeni odnos dve susedne faze kretanja bloka, koje su od posebnog značaja za proces trzanja, i, stoga, nastanak nestabilnosti trenja. Preciznije, vremensko kašnjenje T daje vezu trenja u stenskoj masi sa ,,istorijom'' kretanja bloka pre nego što se blok zaustavio. Direktna posledica ove zavisnosti procesa ,,starenja'' raseda od prethodnog stanja kontaktne zone ogleda se u početnom periodu relativno nepromenjene vrednosti statičkog koeficijenta trenja (čvrstoće na smicanje). Naime, do porasta čvrstoće na smicanje dolazi tek nakon određenog vremenskog intervala od trenutka zaustavljanja bloka. Ovo ponašanje je utvrđeno u laboratorijskim opitima na blokovima granita, koji su klizali preko kataklastičnog, zdrobljenog tektonskog materijala kvarcnog sastava [Marone, 1998]. Opiti su sprovedeni za brzine bloka Vb=0,01-10m/s i pri stalnoj brzini celog sistema (VL=30mm/god.~3x10\"9m/s). Na slici 1-3a jasno se uočava početni period ,,stagniranja'' u porastu čvrstoće na smicanje. Za Vb=Vi=1p,m/s, statička čvrstoća se linearno povećava sa logt. Međutim, za Vb=1m/s, i VL=10-9m/s efektivna brzina regenerisanja raseda ,,kasni'' u vremenskom intervalu odmah nakon dogođenog potresa. Dužina ovog kašnjenja zavisi od krutosti materijala u rasednoj zoni i od veličine rupture koje nastaje usled pojave potresa. Parametri trenja a=0,01; b=0,02; i flc=5mm, u saglasnosti su sa laboratorijskim i terenskim procenama [Marone, 1998]. Pored ovih laboratorijskih rezultata, slično ponašanje je osmatrano i u slučaju realnih raseda, gde se takođe javlja početni period kašnjenja u povećanju čvrstoće na smicanje (<100 dana od poslednjeg zemljotresa), sa malim promenama u seizmičkom stresu (naponu), reda veličine 1MPa, nakon čega dolazi do povećanja napona duž raseda, reda veličine nekoliko MPa (slika 1-3b).\r\n4) uvođenjem vremenskog kašnjenja u izraz za trenje samo se pojačava efekat kompleksne prirode trenja na sam proces kretanja bloka (a i proces kretanja duž raseda in situ). Naime, u rasednoj zoni, usled kretanja, dolazi do pojave zagrevanja i topljenja stenske mase, termalnih pritisaka, stvaranja silika-gela koji deluje kao lubrikant (,,podmazuje rasednu zonu''), kao i efekata granularnog trenja (trenja između zdrobljenih delova stenske mase). U drugom modelu (Madariaga), vremensko kašnjenje T se može smatrati dopunom parametru stanja 6. Naime, vrlo često se u modelima koji uzimaju u obzir Diterih-Ruina zakon trenja, javljaju dva ili više parametara stanja. U ovom slučaju, koristi se samo jedan parametar stanja, ali uz dopunu vremenskim kašnjenjem. Na taj način, pojačavaju se tzv. memorijski efekti kontaktne površine, koja ,,kao da pamti'' prethodne cikluse kretanja, tako da pri svakom novom kretanju, kontaktna površ ima nova, drugačija svojstva;\r\n5) još jedna motivacija za uvođenje vremenskog kašnjenja T proizilazi iz sličnosti BK modela sa relaksacionim oscilatorom. Naime, trzajuće kretanje predstavlja inherentnu karakteristiku relaksacionih oscilatora, čija je glavna odlika postojanje dve vremenske skale - ,,brze'' tokom koje dolazi do oslobađanja akumulirane energije (,,slip'' faza) i ,,spore'' tokom koje se akumulira napon u stenskoj masi (,,stick'' faza). S obzirom na to da je vremensko kašnjenje T inherentno svojstvo relaksacionih oscilatora, u ovom slučaju se uvodi i u BK model, i razmatra se uticaj ovog parametra na dinamiku modela;\r\nSlika 1-3. Vremensko kašnjenje T, osmatrano u laboratorijskim uslovima (a), i in situ (b). Rezultati numeričke simulacije i poređenje sa terenskim podacima [Marone, 1998].\r\n6) u radu [Kijko, i dr., 1993] ukazano je na vremensku zavisnost rudarski izazvanih potresa, u smislu da se parametri stenske mase menjaju tokom vremena, pre i posle seizmičkog događaja. Takva njihova promena uslovljena je rastojanjem od epicentra potresa, veličinom potresa i svojstvima i stanjem stenske mase. S tim u vezi, razmatranje trenja kao funkcije sa usporenim argumentom, omogućava bolje sagledavanje procesa tokom i nakon potresa.\r\nNeophodno je naglasiti da je ovo prvi slučaj uvođenja parametra T u BK model u dosadašnjim istraživanjima. Iako svesni da njegova prava priroda možda i nije egzaktno definisana, imajući u vidu uticaj ovog parametra na dinamiku modela, smatramo da će njegova tačna uloga u realno osmatranim potresima biti određena budućim empirijskim istraživanjima (in situ ili u laboratorijskim uslovima).\r\nRazlozi za uvođenje parametra jačine trenja c u Karlson-Langer model slične su prirode kao i za uvođenje vremenskog kašnjenja T, s tom razlikom da variranjem parametra c, dolazi samo do direktnog smanjenja ili povećanja uticaja trenja (trenje se pojačava ili smanjuje). Ova varijacija jačine trenja u Zemljinoj kori obično je uzrokovana nastankom i svojstvima zdrobljenog materijala duž rasedne zone, koji predstavlja nekonsolidovani materijal bez kohezije, gotovo u potpunosti sastavljen od fino zdrobljenih delova, kao finalnog produkta trenja na kontaktu dva bloka stenske mase [Sibson, 1977]. Postoje dokazi do ovaj materijal ima značajan efekat na stabilnost prirodnih raseda [Scholz, i dr., 1969; Marone, Scholz, 1988] i trenja stenske mase [Byerlee, 1967; Engelder, i dr., 1975; Byerlee, Summers, 1976; Das, dr. 1986], uzrokujući promenu čvrstoće na smicanje duž raseda, u zavisnosti od mineralnog sastava, pre svega prisustva minerala glina u zasićenom stanju [Morrow, i dr., 2000; Behnsen, Faulkner, 2012], zatim u zavisnosti od debljine rasedne zone [Marone, i dr., 1990], kao i od dubine raseda [Mizoguchi, i dr., 2007].\r\nS druge strane, ove promene u jačini trenja mogu biti izazvane i dejstvom nekih spoljašnjih faktora - u ovom slučaju usled napredovanja rudarskih radova. Naime, promenom naponskog stanja, menja se i uticaj trenja duž raseda na početak i sam proces kretanja.\r\nPeriodična perturbacija parametara u Karlson-Langer i Bekerovom sistemu ima za cilj modelovanje spoljašnjih uticaja, u ovom slučaju promena naponskog stanja usled napredovanja rudarskih radova (vibracije rudarskih mašina, efekti miniranja, i dr.). U ovom slučaju, spoljašnje perturbacije napona se uvode implicitno, pretpostavljujući da uzrokuju male oscilacije oko ravnotežnih vrednosti nekih od parametara sistema. U opštem slučaju, ovakvi spoljašnji uticaji mogu se pripisati efektima plime i oseke ili velikih akumulacija [Perfettini i dr., 2001], ili, pak, udaljenih zemljotresa [Gomberg i dr., 1998; Parsons, 2005].\r\nOpšte posmatrano, disertacija se sastoji iz sledećih poglavlja. U prvom delu, u poglavljima ,,Pregled rezultata prethodnih istraživanja'', ,,Mehanizam nastanka rudarski generisanih potresa'', i ,,Modelovanje mehanizma nastanka rudarski generisanih potresa duž raseda'' dat je opšti presek stanja istraživanja u ovoj oblasti, sa prikazom relevantnih podataka o rudarskim potresima i načinu njihovog modelovanja. U okviru poglavlja ,,Nelinearna dinamika i teorija haosa-pregled osnovnih pojmova'' dat je pregled opštih pojmova iz ove oblasti. U poglavlju ,,Metodologija istraživanja'' dat je prikaz korišćenih metoda, dok se poglavlje ,,Rezultati istraživanja'' sastoji iz dva dela: ,,Analiza vremenskih serija rudarski generisanih potresa'' i ,,Analiza dinamike Baridž- Knopof modela''. U poglavlju „Zaključak\" dat je kritički osvrt na dobijene rezultate, sa predlozima daljih istraživanja.\r\n2. PREGLED REZULTATA PRETHODNIH ISTRAZIVANJA\r\nU ovom poglavlju dat je samo opšti pregled rezultata najznačajnijih istraživanja u okviru predmetne naučne oblasti. Detaljan opis rezultata pojedinih istraživanja prikazan je, po potrebi, u okviru svakog poglavlja disertacije.\r\nIstraživanja u oblasti rudarski generisanih potresa uglavnom su zasnovana na proceni hazarda i rizika od njihovog dejstva na sam proces iskopa, zatim napredovanja rudarskih radova, ekonomičnosti eksploatacije, kao i na primenjenu mehanizaciju, objekte rudnika na površini terena i ljudstvo [Gibowicz, Kijko, 1994]. Matematičko modelovanje mehanizma nastanka potresa i njihovih efekata u terenu bili su do sada predmet malobrojnih razmatranja.\r\nPrva veza između nelinearne dinamike i haosa i rudarski izazvanih potresa data je u radu [Morrison, Swan, Scholz, 1993], gde autori, metodom diskretnih elemenata, izučavaju pojavu kretanja duž raseda (a samim tim i potresa) u zavisnosti od rastojanja od čela iskopa. Rezultati njihovog istraživanja upućuju na sledeće zaključke: da raspored potresa (kretanja duž raseda) ima samoslična (fraktalna) svojstva, kao i da je deterministički model u slučaju potresa prihvatljiv samo do izvesne granice, odnosno da je poboljšanje modela moguće ukoliko pođemo od pretpostavke da je nastanak potresa ,,haotični'' sistem.\r\nKortas [2005] istraživao je seizmograme 5 jačih potresa zabeleženih u Poljskoj i Južnoj Africi. Dva su registrovana u basenu bakra Legnica-Glogow (u rudnicima Polkovice-Sieroszovice i Rudna), dva u basenu uglja u Gornjoj Sileziji (u rudnicima Katovice i Vujek), dok peti seizmogram predstavlja registrovane oscilacije tla za vreme potresa u jednom rudniku zlata u Južnoj Africi. Istraživanjem je obuhvaćen raspored registrovane energije, vremena registrovanja, koordinata epicentara i hipocentara potresa. Koristeći Grasberger-Prokačia algoritam za određivanje veličine fraktalne dimenzije, Kortas je ukazao na to da distribucija vremena pojave, energije i koordinata epicentara ima svojstvo multifraktala (0,69-1,60 za rudnike bakra, 0,62-2,56 za rudnike uglja, i 0,87-2,07 za rudnik zlata u Južnoj Africi). Takođe, Kortas je pristupio izračunavanju vrednosti maksimalnog Ljapunovljevog eksponenta, pri čemu je u 4 od proučavanih 20 slučajeva (pet rudnika sa četiri različite distribucije registrovanih parametara) maksimalni Ljapunovljev eksponent imao pozitivnu vrednost, ukazujući na visokodimenzionalno determinističko haotično ponašanje sistema.\r\nS druge strane, Huang, Yin i Dai [2010] istraživali su model dva bloka sa oprugom, u cilju simulacije rudarski generisanih potresa, i pri tome utvrdili determinističko haotično ponašanje, koje ukazuje na to da je dinamika modela rudarski generisanih potresa deterministička, nasuprot često zastupanom stavu o stohastičkoj dinamici modela potresa.\r\nNeophodno je naglasiti da je primena ovih metoda vrlo zastupljena u ispitivanju tektonskih, prirodnih zemljotresa. Imajući u vidu činjenicu da je sam proces kretanja duž raseda, ukoliko zanemarimo razlike u načinu pobude, gotovo isti kod tektonskih i rudarski generisanih potresa, postojeća saznanja iz oblasti seizmologije se analogno primenjuju i u oblasti rudarske seizmologije1 [Gibowicz, Kijko, 1994], poput Gutenberg-Rihter zakona, fokalnih mehanizama i fenomenoloških modela. S obzirom na to da ovaj deo istraživanja obuhvata analizu fenomenoloških modela, u daljem tekstu daćemo pregled rezultata dosadašnjih istraživanja u ovoj oblasti.\r\nEmpirijska osnova za kreiranje fenomenološkog BK modela postavljena je u radu [Brace, Byerlee, 1966], pri čemu su autori vršili ispitivanja čvrstoće na smicanje poremećenog i neporemećenog uzorka granita u triaksijalnom aparatu. Rezultati njihovih ispitivanja ukazali su na to da u neporemećenom uzorku, kada vertikalna sila prevaziđe čvrstoću na smicanje stenske mase i dođe do pojave loma, dalje kretanje jednog dela uzorka u odnosu na drugi postaje trzajuće2. Naime, neposredno nakon loma, duž kontaktne površi dva dela uzorka formira se jedna zona zdrobljene stenske mase, sa velikom vrednošću koeficijenta trenja, odnosno čvrstoće na smicanje. Daljim delovanjem sile u triaksijalnom aparatu, dolazi do povećanja čvrstoće na smicanje, pri čemu nema kretanja duž kontaktne površi dva dela uzorka. Onog momenta kada sila u triaksijalnom aparatu ,,prevaziđe'' veličinu statičkog koeficijenta trenja, dolazi do ponovnog kretanja. Imajući u vidu ove rezultate, u radu [Brace, Byerlee, 1966] predloženo je da bi ovakvo kretanje moglo da opiše nastanak zemljotresa u Zemljinoj kori. Nedugo zatim, oslanjajući se na ove rezultate, Baridž i Knopof [Burridge,\r\nKnopoff, 1967] konstruisali su model bloka sa oprugom, koji je oponašao trzajuće kretanje, pri čemu su pretpostavili da je zakon trenja duž kontaktne površi funkcija samo brzine. Detaljan opis BK modela dat je u poglavlju ,,Mehanizam nastanka rasednih potresa''. Do sličnih rezultata došli su i Karslon i Langer [Carlson, Langer, 1989]. Međutim, naredna istraživanja su pokazala da je trenje duž kontakta dva bloka stenske mase isuviše kompleksno da bi se moglo opisati jednostavnim zakonom trenja zavisnim samo od jedne promenljive (brzine). Shodno tome, Diterih [Dieterich, 1979] je izveo niz laboratorijskih opita, čije je rezultate opisao čuvenim Diterih zakonom trenja, gde je trenje zavisno ne samo od brzine, već i od stanja površi duž koje se blok kreće. Ovaj zakon trenja je danas najčešće poznat kao Diterih-Ruina zakon trenja, zbog značajanog doprinosa A. Ruine [Ruina, 1983] daljem matematičkom razvijanju predloženog zakona.\r\nSa razvojem brzih računara složenije numeričke kalkulacije postale su moguće, što je i otvorilo put za primenu znanja iz oblasti nelinearne dinamike u istraživanju zemljotresa. Tako još Karlson i Langer [Carslon, Langer, 1989] ukazuju na mogućnost pojave ,,haotičnog'' ponašanja u BK modelu. Potom, Vieira [Vieira, 1999] dokazuje prisustvo determinističkog haosa u modelu sa jednim i tri bloka, pretpostavljajući zakon trenja zavisan od brzine. Beker [Becker, 2000] proučavajući model jednog bloka sa oprugom, i zakonom trenja zavisnim od brzine i stanja, dolazi da zaključka da varirajući krutost opruge dolazi do pojave determinističkog haosa i to preko niza bifurkacija sa udvajanjem perioda (Fajgenbaumov put u haos). Erikson i dr. [Erickson, i dr., 2008], na osnovu analize Madariaga sistema, koji matematički opisuje kretanje BK modela sa jednim blokom, takođe pokazuju da za određene vrednosti parametara dolazi do pojave deterministički haotičnog ponašanja. Takođe, isti autori [Erickson, i dr., 2011] dalje pokazuju, da, u slučaju diskretnog BK modela, sistem ispoljava i periodično i ,,haotično'' ponašanje, pri čemu je prelaz ka ,,haotičnom'' ponašanju uslovljen brojem posmatranih blokova. Iz diskretnog modela autori dalje izvode nelinearnu talasnu jednačinu, čijim rešavanjem se ukazuje na postojanje prostorno-vremenskog haotičnog ponašanja modela. Rezultati istraživanja upućuju na zaključak da je vrednost kritičnog parametra potrebna za pojavu haosa u kontinualnom modelu mnogo manja od vrednosti u slučaju samo jednog bloka. Oni takođe dolaze do zaključka da je usvojeni Diterih- Ruina zakon trenja zavisan od posmatrane razmere.\r\nS druge strane, haotično ponašanje osmatrano u fenomenološkim modelima, ispitivano je i kod vremenske distribucije registrovanih potresa. Tako su Beltrami i Marešal [Beltrami, Mareschal, 1993] pokušali da rekonstrišu strani atraktor za vremensku seriju registrovanih potresa u seizmičkoj oblasti Parkfild u Kaliforniji, u periodu 1969-1987.g. Pri tome, došli su do zaključka da je serija potresa ili potpuno slučajna, ili ima strani atraktor (,,haotična''), sa dimenzijom atraktora većom od 12. Tivari i dr. [Tiwari, i dr., 2004] primenili su postupak nelinearnog predviđanja u rekonstruisanom faznom prostoru pri proučavanju frekvencije pojave potresa u oblasti Centralnih Himalaja, i pri tome su utvrdili nisku pozitivnu korelaciju između predviđenih i realno osmatranih podataka, ukazujući na zajedničku stohastičku i ,,haotičnu'' prirodu. Rezultati Iljopulosa i dr. [Iliopulous, i dr. 2008] ukazuju na prostorno-vremensku ,,haotičnu'' prirodu seizmičkog procesa, i na postojanje globalnog ,,zemljotresnog'' stranog atraktora, kao osnovnog fizičkog procesa seizmogeneze u oblasti Grčke. De Santis i dr. [De Santis, i dr., 2010] pokazali su da se niz zabeleženih prethodnih udara, koji su kulminirali potresom magnitude Mw = 6.3 u oblasti Akvila u centralnoj Italiji 2009.g., razvijao kao ,,haotični'' proces. Pri tome, primenili su metodu zasnovanu na ubrzanom oslobađanju napona u vremenu3, uz nelinearni pristup u rekonstruisanom faznom prostoru.\r\nPojedini autori imaju različite stavove u pogledu uzroka pojave determinističkog haosa u zabeleženim vremenskim serijama. Tako jedna grupa autora smatra da je deterministički haotična priroda posledica samog mehanizma nastanka potresa [Vieira, 1992; Erickson, i dr., 2008; Carlson, Langer, 1989], dok drugi smatraju da se ,,haos'' javlja usled propagacije seizmičkih talasa kroz heterogene stenske mase [Abdullaev, Zaslavskii, 1991; McCloskey i dr., 1991; Urquizu, Correig, 1998].\r\nNeophodno je naglasiti da koncept potresa kao deterministički haotičnih pojava nije podržan od strane pojedinih autora. Tako pojedini autori smatraju da seizmičnost predstavlja čisto slučajan (nasumičan) proces [Rundle, i dr., 1989; Huang, Turcotte, 1990; Kagan, 1997; Wyss, Toya, 2000], dok druga grupa autora smatra da su zemljotresi u stvari predstavnici klase samo-organizovanih kritičnih pojava [Bak, i dr., 1988; Malamud, Turcotte, 2000; Lise, Paczuski, 2001; Caruso, i dr., 2007].\r\nU pogledu analize zabeleženih oscilovanja tla, metode nelinearne analize i analize surogat podataka nisu do sada primenjivane u oblasti rudarskog inženjerstva. S druge strane, pojedine metode su primenjivane u istraživanju oscilovanja tla tokom tektonskih zemljotresa. Tako su Pavlos, i dr. [1994], analizirajući vremenske serije zemljotresa na području Japana, ustanovili da mehanizam nastanka potresa, koji je opisan registrovanom vremenskom serijom, predstavlja niskodimenzionalni haotični proces. Takođe, u radu [Yang et al., 2012] izučavane su vremenske serije zabeleženog oscilovanja tla tokom Či-Či zemljotresa, magnitude 7,6, koji je 1999.g. pogodio centralni Tajvan. Pri tome, analiza registrovanih potresa na 24 stanice ukazala je na ,,haotičnu'' i fraktalnu prirodu oscilovanja tla tokom potresa.\r\n3. MEHANIZAM NASTANKA RUDARSKI GENERISANIH POTRESA\r\nPotresi nastali rudarskom aktivnošću pri iskopu podzemnih prostorija za potrebe eksploatacije mineralnih sirovina, ili u toku same eksploatacije, spadaju u grupu zemljotresa tehnogenog porekla. Pri tome, razlikuju se dve osnovne grupe potresa: indukovani i pobuđeni, u zavisnosti od seizmičnosti datog područja. Ukoliko se potresi usled rudarske aktivnosti1 jave u seizmički aktivnim područjima, onda je reč o pobuđenim zemljotresima. U suprotnom, u slučaju pojave potresa u aseizmičkim područjima, zemljotresi se nazivaju indukovanim.\r\nSeizmičnost izazvana dubokim podzemnim radovima u rudarstvu (eksploatacijom mineralnih sirovina, izgradnjom velikih podzemnih prostorija) predstavlja najznačajniji geološki hazard, koji direktno utiče kako na sam proces eksploatacije, tako i na bezbednost ljudi i opreme.\r\nIzradom podzemnih prostorija za potrebe eksploatacije mineralnih sirovina, kao i tokom samog procesa otkopavanja rude, dolazi do poremećaja primarnog naponskog stanja, pri čemu su novonostale promene veće i složenije, ukoliko su inženjerskogeološki uslovi u stenskoj masi pre izrade podzemne prostorije bili složeni, a takođe i ukoliko proces otkopavanja, kao i veličina i oblik prostorija nisu prilagođeni prirodnim uslovima koji vladaju u tom delu stenske mase [Cabala, i dr., 2004]. Ova preraspodela napona, u pojedinim slučajevima, može dovesti do pojave potresa, pri čemu potresi najveće magnitude (M>4) po pravilu bivaju izazvani preraspodelom tektonskih napona u neposrednoj blizini otkopa. Prethodna istraživanja [Smith, i dr., 1974] ukazuju na to da rudarska aktivnost u Zemljinoj kori, u opštem slučaju, dovodi do smanjenja veličine litostatičkog pritiska, putem rasterećenja stenske mase. Istovremeno, horizontalni naponi su takođe izmenjeni, kao posledica rasterećenja i Puasonovog efekta [Wong, 1985]. U kompresionom tektonskom polju, koje preovladava u Zemljinoj kori, efekat rudarske aktivnosti izazvaće kretanja duž postojećih mehanički oslabljenih zona (raseda) pre pojave loma, i to smanjenjem vertikalno orijentisanog maksimalnog glavnog napona, i povećanjem smičućeg napona duž raseda. Kao rezultat toga, može doći do stvaranja seizmički aktivne rupture, u zavisnosti od veličine koeficijenta trenja duž pukotina, i njihove orijentacije u odnosu na pravac maksimalnog glavnog napona.\r\nU opštem slučaju, rudarski izazvani potresi javljaju se u područjima sa srednjim do visokim horizontalnim kompresionim tektonskim poljem, koji se istovremeno odlikuju srednjom do visokom tektonskom seizmičnošću. Najjači potresi među njima, za koje su dostupni fokalni mehanizmi, u opštem slučaju nastaju kao rezultat reverznog ili transkurentnog rasedanja. U znatno manjem broju slučajeva, i normalni tip rasedanja može dovesti do pojave potresa (slika 3-1). Ovakav tip deformacije predstavlja posledicu povećanja napona smicanja duž postojećih oslabljenih zona (raseda), usled rasterećenja stenske mase. U slučaju kada litostatički pritisak određuje naponsko stanje, tektonski naponi mogu delovati kao ,,okidač'' za pojavu rudarski izazvanih potresa.\r\nSlika 3-1. Shematizovani prikaz osnovnih mehanizama nastanka rudarski izazvanih potresa usled kretanja duž raseda (promena napona): a) kretanje duž normalnog raseda; b) kretanje duž reversnog raseda; c) kretanje duž transkurentnog raseda [Hasegawa, i dr., 1989].\r\nNeophodno je napomenuti da je lom smicanjem najčešći uzrok nastanka potresa i u aseizmičkim oblastima, odnosno područjima u kojima prethodno nisu zabeleženi tektonski zemljotresi. Tada, usled redistribucije napona u netaknutoj stenskoj masi, u prvoj fazi dolazi do formiranja tenzionih pukotina, s obzirom na to da je čvrstoća na zatezanje najmanja kod stenskih masa. U sledećoj fazi, preraspodela napona prevazilazi i čvrstoću na smicanje, kada dolazi do pojave pukotina smicanja. Na kraju, daljim napredovanjem otkopavanja, dolazi do kretanja duž formirane pukotine smicanja, pri čemu akumulirana potencijalna energija biva delom transformisana u kinetičku, a delom se oslobađa putem toplote i seizmičkih talasa, kada se stvaraju uslovi za pojavu potresa.\r\nPromena naponskog stanja, pored toga što direktno može da dovede do pojave potresa, često uzrokuje pojavu deformacija stenske mase, kada takođe dolazi do slabijih potresa (slika 3-2). Naime, kao posledica deformacija stenske mase može doći do gorskih udara2, razaranja zaštitnih stubova3, iznenadnog klizanja kvazi-viskoznog rudonosnog sloja na velikim dubinama4 (pri čemu ovaj termin obuhvata i potrese nepoznatog porekla), kao i do naglog izbacivanja ugljenog proslojka u otkopani prostor5. Pored toga, do potresa izazvanih deformacijama stenske mase može doći i primenom miniranja na čelu otkopa, kada usled posledica razaranja stenskog materijala i brze konvergencije stvorenog prostora dolazi do formiranja zdrobljenih zona klizanja6, duž kojih može doći do kretanja stenske mase i pojave potresa.\r\nSlika 3-2. Shematizovani prikaz osnovnih mehanizama nastanka rudarski izazvanih potresa usled deformacije stenske mase: a) prolom krovine podzemnih prostorija; b) lom zaštitnih stubova; c) formiranje tenzionih pukotina u krovini. Strelice označavaju pravac dejstva napona (stresa) [Hasegawa, i dr., 1989].\r\nOvakva podela potresa izazvanih izradom horizontalnih podzemnih prostorija predstavlja opštu podelu, i moguće je izdvojiti nekoliko podtipova. Šest modela rudarski indukovane seizmičnosti, originalno razvijenih na osnovu zabeleženih potresa u rudnicima u Kanadi [Hasegava, i dr. 1989], prihvaćeni su kao opšti tipovi mehanizma nastanka potresa izazvanih podzemnim radovima:\r\n- potresi izazvani zarušavanjem krovine tokom procesa otkopavanja;\r\n- potresi nastali usled loma zaštitnih stubova;\r\n- potresi nastali usled formiranja tenzionih pukotina u krovini podzemnih prostorija;\r\n- potresi izazvani kretanjem duž normalnog (gravitacionog) raseda;\r\n- potresi nastali pomeranjem duž reversnog raseda;\r\n- potresi izazvani kretanjem duž transkurentnog raseda.\r\nPored ove klasifikacije mehanizma nastanka rudarski generisanih potresa, u radu [Ortlepp, 1992] predloženo je pet osnovnih tipova potresa prema mehanizmu nastanka (Tabela 3-1), u zavisnosti od tipičnog kretanja stenske mase, pravca kretanja koji se najpre registruje instrumentalno, kao i približnog raspona magnituda potresa.\r\nPredloženi mehanizmi odgovaraju lomovima u stenskoj masi velikih razmera. U opštem slučaju, međutim, više od 90% seizmičke aktivnosti predstavlja manje lomove u stenskoj masi od onih prethodno opisanih od strane Hasegave i Ortlepa. Takvi mehanizmi mikroseizmičkih događaja (mikropotresa) u rudnicima obuhvataju sve prethodne mehanizme predložene u radu [Ortlepp, 1992], ali takođe i sledeće mehanizme:\r\n- krti lom intaktne stenske mase [Lynch, i dr., 2005];\r\n- srastanje pukotina u stenskoj masi [Trifu, Urbancic, 1996];\r\n- koncentraciju napona na krajevima otkopa [Andrieux, Simser, 2001];\r\n- lom zaštitnih stubova drobljenjem i smicanjem [Hedley, 1992], i\r\n- smicanje duž kontakta litološki različitih vrsta stenskih masa [Mollison i dr., 2001].\r\nTabela 3-1. Osnovni mehanizmi nastanka razornih potresa izazvanih rudarskom aktivnošću [Ortlepp, 1992].\r\n*negativna vrednost magnitude označava događaje vrlo male razmere (s obzirom na to da se magnitude izražavaju u logaritamskoj razmeri u odnosu na amplitudu, desetostruko smanjenje amplitude odgovara smanjenju magnitude od 1; tako da ukoliko amplituda kretanja tla od 20mm odgovara magnitudi 2, onda će 1000 puta manja amplituda (0,02mm) odgovarati magnitudi od -1; intuitivno je jasno da se potresi ove jačine instrumentalno ne mogu registrovati).\r\nU radu [Hudyma, 2004] dat je ilustrativni primer mogućih mehanizama nastanka rudarski generisanih potresa, sa priloženim fotografijama realnih događaja in situ (slika 3-3).\r\nPrethodna istraživanja su pokazala da su potresi u rudnicima izazvani reaktivacijom kretanja duž raseda jedna od najčešćih pojava rudarski indukovane seizmičnosti [Donnelly, Reddish, 1994; Donnelly, 2006; Cui, i dr. 2004; Li, i dr. 2004]. Štaviše, ovi potresi uglavnom predstavljaju najveće događaje sa najtežim posledicama. Najjači potresi često dostižu magnitudu M=5, u izuzetnim slučajevima i do M=6. Udari ovakvih magnituda sa dubinom žarišta od 1-2km mogu izazvati potrese intenziteta VI-VIII stepeni po Modifikovanoj Merkalijevoj skali [Zembaty, 2011].\r\nSlika 3-3. Opšti primeri loma u stenskoj masi koji mogu uzrokovati potrese: (a) kretanje duž raseda; (b i f) promena napona koja uzrokuje deformaciju stenske mase u blizini iskopa; (c) lom na krajevima čela otkopa; (d) razlika u litološkom sastavu stenske mase uzrokuje gorski udar; (e) lom zaštitnih stubova [Hudyma, 2007; prilagođeno].\r\n3.1. RASEDNI POTRESI12\r\nKao što je već pomenuto, ovakav tip rudarski izazvanih potresa je najčešći i sa najvećim magnitudama. U ovom slučaju, gotovo da nema razlike između rudarski generisanih i tektonskih zemljotresa u pogledu mehanizma nastanka - u oba slučaja do pojave potresa dolazi usled loma smicanjem u krtoj stenskoj masi [McGarr, 1971; Spottiswoode, McGarr, 1975; Wong, 1993). U radu [Brace, Byerlee; 1966] prvi put je predloženo ,,trzajuće kretanje''13 duž postojećih raseda kao glavni uzrok nastanka zemljotresa, zahvaljujući rezultatima laboratorijskih opita čvrstoće na smicanje na poremećenom i neporemećenom uzorku granita u triaksijalnom aparatu. U radu [McKenzie, 1969] tvrdi se da do kretanja duž postojećih raseda dolazi pri veoma niskim vrednostima smičućih napona, mnogo pre loma stenske mase i formiranja novih prslina i pukotina. U radu [Rayleigh, i dr., 1972] predloženo je da do seizmičkog kretanja duž postojećih oslabljenih zona dolazi pre loma intaktne stenske mase, ukoliko je ta oslabljena zona orijentisana pod uglom od 10-50° u odnosu na pravac maksimalnog glavnog napona (slika 3-4). Ovaj zaključak je proistekao iz razmatranja pojave indukovane seizmičnosti usled injektiranja fluida na velikim dubinama na Rejndžli naftnom polju u Koloradu.\r\nSlika 3-4. Blok dijagram reverznog raseda u polju kompresionog napona, gde su maksimalni (o^) i srednji glavni napon (o2) horizontalni. Minimalni glavni napon (o3), usled dejstva litostatičkog pritiska je vertikalan. Do kretanja duž raseda dolazi pre loma intaktne stenske mase, ukoliko je rased pod uglom 10-500 u odnosu na pravac maksimalnog glavnog napona [Wong, 1993; prilagođeno].\r\n3.2. PROMENA NAPONA I POMERANJA DUŽ RASEDNE POVRŠI SA NAPREDOVANJEM RUDARSKIH RADOVA\r\nRezultati istraživanja potresa koji se dogodio 30. novembra 2004.g u rudniku uglja u oblasti Jining u Kini [Li i dr., 2008], ukazuju na značajan uticaj rudarskih radova na promenu napona duž rasedne površi, u zavisnosti od toga da li rudarski radovi napreduju sa podinskog krila (slika 3-5a) ili sa povlatnog krila raseda (slika 3-5b).\r\nSlika 3-5. Raspored sila duž rasedne površi pri napredovanju rudarskih radova: (a) sa podinskog krila raseda; (b) sa povlatnog krila raseda. [Li i dr., 2008; prilagođeno]\r\nKada se iskop vrši sa podinskog krila raseda (slika 3-5a), normalna sila duž rasedne površi data je izrazom: N=Tcos6-Rsin6, a sila smicanja duž rasedne površi izrazom: F=Rcos6+Tsin6. U ovom slučaju, uslov stabilnosti (ravnoteže) glasi:\r\nodnosno\r\ngde je T horizontalna sila, R vertikalna sila, P ugao unutrašnjeg trenja, a 0 padni ugao padne prave rasedne površi.\r\nKada iskop napreduje sa povlatnog krila, normalna sila duž rasedne površi data je izrazom N=Tcos0+Rsin0, a smičuća sila duž rasedne površi: F=Rcos0-Tsin0. U ovom slučaju, uslov stabilnosti (ravnoteže) glasi:\r\nSa napredovanjem rudarskih radova, dolazi do konvergencije podzemne prostorije, i krovina se polako spušta prema podu podzemne prostorije. Kako radovi napreduju, smanjuje se oslonac krovine u vidu ugljenog sloja, i pri tome dolazi do povećanja vertikalne sile R. Ukoliko zanemarimo promene u horizontalnoj sili T, onda u prvom slučaju (slika 3-5a), normalna sila se smanjuje, dok se smičuća sila povećava. U ovim uslovima, može doći do pojave kretanja duž raseda, a samim tim i do pojave potresa. S druge strane, kada iskop napreduje iz povlatnog krila, tada do pojave nestabilnosti dolazi takođe usled povećanja vertikalne sile R, a kada je čelo iskopa u blizini rasedne površi, tada dolazi do koncentracije napona duž raseda. Takođe, tada se stvaraju uslovi i za pojavu potresa.\r\nRezultati dobijeni numeričkom analizom geološkog modela bušotine C8-9, u rudniku uglja u oblasti Jining, ukazuju na različiti karakter promene napona, u zavisnosti od rastojanja čela iskopa od rasedne površi. Naime, u slučaju kada rudarski radovi napreduju sa podinskog krila raseda, promena odnosa smučućeg i normalnog napona u zavisnosti od rastojanja iskopa od rasedne površi prikazana je na slici 3-6.\r\nSlika 3-6. Promena odnosa smičućeg i normalnog napona u zavisnosti od rastojanja čela iskopa od rasedne površi. Napredovanje radova je sa podinskog krila. [Li i dr., 2008; prilagođeno]\r\nKao što se sa slike 3-6 može videti, kada je normalni napon oko 30MPa, a smičući napon 2MPa, čelo iskopa je na velikoj udaljenosti od raseda (40-80m). Normalni napon počinje da opada, a smičući napon polako raste, kada je čelo iskopa na rastojanju 30m od rasedne površi. Potom, smičući napon se naglo povećava kada je čelo iskopa oko 20m od raseda, pri čemu je veličina normalnog napona jednaka polovini početne vrednosti. Međutim, kada je čelo iskopa na udaljenosti oko 10m, smičući napon se takođe smanjuje, ali tada odnos smičućeg i normalnog napona pokazuje vrlo visoke vrednosti. Tada dolazi i do kretanja duž raseda, odnosno do pojave rudarski izazvanih potresa.\r\nS druge strane, kada radovi napreduju sa povlatnog krila, dijagram zavisnosti odnosa smičućeg i normalnog napona u funkciji rastojanja od rasedne površi izgleda drugačije i prikazan je na slici 3-7.\r\nSlika 3-7. Promena odnosa smičućeg i normalnog napona u zavisnosti od rastojanja čela iskopa od rasedne površi. Napredovanje radova je sa povlatnog krila. [Li i dr., 2008; prilagođeno]\r\nU ovom slučaju, normalni napon se postepeno povećava, dok se smičući napon najpre povećava, a potom počinje da opada. Kada je udaljenost čela iskopa od rasedne površi preko 30m, dolazi do povećanja i normalnog i smičućeg napona, ali njihov odnos ostaje gotovo konstantan i iznosi oko 0,09.\r\nNumerička ispitivanja obuhvatila su takođe i procenu veličine pomeranja duž raseda, u zavisnosti od udaljenosti čela iskopa od rasedne površi. Na slici 3-8 prikazano je\r\npoređenje veličine pomeranja duž rasedne površi, za slučajeve kada radovi napreduju sa podinskog i povlatnog krila raseda. Dve krive se gotovo podudaraju na rastojanju većem od 30m od čela iskopa do rasedne površi. Sa napredovanjem rudarskih radova, pomeranje duž raseda se u oba slučaja postepeno povećava. Međutim, sa daljim približavanjem iskopa rasednoj površi dolazi do naglog povećanja pomeranja za slučaj kada radovi napreduju iz podinskog krila raseda.\r\nSlika 3-8. Poređenje pomeranja duž rasedne površi u zavisnosti od rastojanja čela iskopa, kada radovi napreduju iz podinskog i povlatnog krila raseda. [Li i dr., 2008; prilagođeno]\r\n3.3. VELIČINA ZONE UTICAJA OKO PODZEMNE PROSTORIJE U KOJOJ MOŽE DOĆI DO POJAVE POTRESA\r\nPromene naponskog stanja obično se razmatraju za slučaj ravnomernog (hidrostatičkog) opterećenja oko nepodgrađene podzemne prostorije, imajući u vidu efekat podgrade jedino u vidu smanjenja konvergencije profila, sprečavanje deformacija stenske mase i pojave podzemnih pritisaka na podgradu. Razmatranje promena naponskog stanja u stenskoj masi, izazvanih iskopom podzemnih prostorija, vrši se svođenjem posmatranog slučaja na neki od prihvaćenih fizičkih modela: elastični, plastični, nevezani, reološki, kao i njihove kombinacije. Međutim, stenska masa u terenu je heterogena i anizotropna sredina, tako da već postojeća rešenja iz teorijske mehanike ne odgovaraju uvek realno osmatranim podacima. U najvećem broju slučajeva, pri izradi horizontalnih podzemnih prostorija na većim dubinama, stenska masa pokazuje krto ponašanje, što uslovljava stvaranje tri zone oko izvedene podzemne prostorije: ispucale, plastične i elastične [Jovanović, 1994]. U praksi, raspored napona oko podzemne prostorije predstavlja kombinaciju ova tri slučaja. U zoni razaranja, deformacije su takvog obima da dolazi do stvaranja pukotina različitog obima i orijentacije, stvarajući oko prostorije jedan nestabilan prstenasti omotač, sklon obrušavanju. U ovoj zoni stenska masa je izdeljena mnogobrojnim pukotinama, tako da je moguće pretpostaviti da je kohezija jednaka nuli. Stenska masa se u ovom slučaju nalazi u stanju labilne ravnoteže, koja se održava samo zahvaljujući uglu unutrašnjeg trenja. Mala zaostala naprezanja koja postoje na konturi zone razaranja se sa približavanjem plastičnoj zoni povećavaju, da bi na granici sa plastičnom zonom tangencijalni normalni napon dobio vrednost koja odgovara vrednosti čvrstoće stenske mase. U zoni plastičnih deformacija naprezanje se i dalje povećava sve do granice sa elastičnom zonom, kada tangencijalni napon dostiže svoj maksimum, a radijalni sa nešto većom brzinom teži da se što pre približi naponu u neporemećenom stenskom masivu. U elastičnoj zoni, tangencijalni napon opada i asimptotski se približava prvobitnoj vrednosti napona, dok radijalni nastavlja da se povećava i asimptotski teži ka početnoj vrednosti napona, koja odgovara naponu kod neporemećenog masiva (slika 3- 9). Neophodno je naglasiti da, iako se u praksi podzemne prostorije izvode sa kružnim, eliptičnim i pravougaonim oblikom poprečnog preseka, rešenja za elasto-plastične uslove data su samo za kružni oblik poprečnog preseka.\r\nVeličina zone uticaja promene naponskog stanja oko horizontalne podzemne prostorije kružnog poprečnog preseka (poluprečnik R), u uslovima ravnomernog (hidrostatičkog) opterećenja za slučaj elasto-plastičnog ponašanja stenske mase, može se dati sledećim obrascem [Jovanović, 1994]:\r\ngde je:\r\nrk - poluprečnik krto-razorene zone, rp - poluprečnik uticaja plastične zone, rk - poluprečnik uticaja elastične zone, op- čvrstoća na pritisak stenske mase, p - ugao unutrašnjeg trenja, Y - zapreminska težina stenske mase,\r\nH - dubina podzemne prostorije u odnosu na površinu terena.\r\nPrema tome, na ovom maksimalnom rastojanju od podzemne prostorije (čela iskopa) može se o čekivati i pojava potresa, ukoliko na toj udaljenosti postoje rasedne zone duž kojih može doći do reaktiviranja kretanja.\r\nSlika 3-9. Shema elasto-plastičnog modela. 1 - zona razorenog materijala, 2 - plastična zona, 3 - elastična zona [Jovanović, 1994; prilagođeno].\r\nLegenda: or - radijalni napon; ot - tangencijalni napon; Rp - jednoaksijalna čvrstoća na pritisak stenske mase; rk - poluprečnik zone razorenog materijala; rp - poluprečnik plastične zone; p - podzemni pritisci; a - poluprečnik podzemne prostorije; u - deformacija konture, koja se javlja kao posledica uravnoteženja naponskog stanja oko podzemne prostorije.\r\n3.4. REGISTROVANE POJAVE RASEDNIH POTRESA U SVETU\r\nZa najveći broj jačih registrovanih potresa, izazvanih rudarskom aktivnošću, mehanizam nastanka može se objasniti reaktiviranjem kretanja duž raseda, ili lomom na smicanje intaktne stenske mase. U daljem tekstu dat je pregled značajnijih rudarskih oblasti u svetu, u kojima su registrovani potresi sa prikazom njihovog mehanizma nastanka.\r\nU Sjedinjenim Američkim državama najbolje dokumentovani slučajevi rasednih potresa registrovani su u rudnicima olova, srebra i cinka u oblasti Kor d' Alen14 u severnom delu države Ajdaho, kao i u rudnicima uglja u oblastima Vosač15 i Buk Klifs u centralnom delu države Juta (slika 3-10). Obe ove oblasti nalaze se na zapadnoj obali Severne Amerike, gde dominira nekoliko režima tektonskih napona [Zoback, Zoback, 1989], pre svega zbog blizine granice tektonskih ploča duž pacifičke obale i visokog stepena seizmičnosti. Fokalni mehanizmi za 106 rudarski generisanih potresa u rudniku Laki Frajdej16 u Kor d' Alen oblasti pokazuju predominatno transkurentno i reversno gravitaciono kretanje duž raseda sa pružanjem sever-jug (slika 3-10). Ovi mehanizmi nastanka u skladu su sa kompresionim naponskim poljem severozapadnog Pacifika. Merenja na jezgru takođe su pokazala da je maksimalni glavni napon u pravcu severozapad-jugoistok oko 1,3 puta veći od litostatičkog ili minimalnog glavnog napona (Lourence, i dr., 1993). Takođe, fokalni mehanizmi rudarski generisanih potresa se u ovom slučaju poklapaju sa mehanizmima tektonskih zemljotresa, što ukazuje na zajedničko, tektonsko poreklo. S druge strane, na Vosač platou, fokalni mehanizmi rudarski generisanih potresa, maksimalne magnitude M=3,5, pokazuju reversno rasedanje, pružanja zapad do severozapad (Wong, i dr., 1989). U ovom slučaju, Vosač Plato se nalazi u blizini intraplaninskog seizmičkog pojasa, u kojem su fokalni mehanizmi takođe pokazali reversno rasedanje. Prema tome, i tektonski zemljotresi i rudarski generisani potresi javljaju se u naponskom polju koje se karakteriše kompresionim naponima u pravcu severoistok-jugozapad.\r\nSlika 3-10. Tipični fokalni mehanizmi rudarski izazvanih potresa u državi Ajdaho i Juta. Strelice orijentisane ka fokalnim mehanizmima označavaju horizontalne projekcije maksimalnog glavnog napona, a strelice orijentisane od fokalnih mehanizama - projekcije minimalnih glavnih napona. Mehanizmi 1-3 preuzeti su iz rada [Lourence, i dr. 1993], a 4-5 iz rada [Wong, i dr., 1989]. Krupne strelice pokazuju pravac regionalnog maksimalnog napona. [Wong, 1993; prilagođeno]\r\nU istočnoj Kanadi, u Sadberi basenu, u dubokim rudnicima metaličnih mineralnih sirovina, registrovan je veći broj rudarski izazvanih potresa. Ovaj deo istočne Kanade takođe je seizmički aktivan, pri čemu fokalni mehanizmi registrovanih zemljotresa pokazuju reversno, a ponekad i transkurentno rasedanje u kompresionom tektonskom naponskom polju [Adams, Bell, 1991]. Maksimalni glavni napon u istočnoj Kanadi orijentisan je u pravcu severoistok do istok, mada može doći i do velikih lokalnih varijacija. Ovi horizontalni kompresioni naponi su vrlo visoki, čak i do dva puta veći od litostatičkog pritiska. U radu [Wetmiller, i dr., 1990] određeni su fokalni mehanizmi za pet većih rudarski izazvanih potresa u Sadberi basenu u periodu 1984-1987.g. (Slika 3- 11). Mehanizmi dva potresa u Straćona rudniku, sa magnitudnim momentom mb=3,1 i 3,0, pokazuju transkurentno i reversno rasedanje, sa pružanjem severozapad-jugoistok i istok-zapad. Fokalni mehanizmi tri potresa registrovana u Krejton rudniku, sa magnitudnim momentom mb= 4,0; 3.3 i 3,7, pokazuju reversno rasedanje. Prema tome, četiri od pet ovih većih registrovanih potresa su jako slični tektonskim zemljotresima, pošto nastaju kao rezultat reversnog rasedanja u kompresionom tektonskom polju istočne Kanade (slika 3-11).\r\nSlika 3-11. Tipični fokalni mehanizmi rasednih potresa u Sadberi basenu u istočnoj Kanadi (Wetmiller et al., 1990). [Wong, 1993; prilagođeno]\r\nZa rudarski generisane potrese u Poljskoj, u radu [Gibowicz, 1984] određeni su fokalni mehanizmi za četiri velika potresa koji su registrovani u tri različite rudarske oblasti u Poljskoj. Dva potresa, prvi, magnitude ML=4,5, koji se desio 24. juna 1977.g. između dva rudnika bakra u oblasti Lubin, i drugi, magnitude ML=4,3, koji se desio 30. septembra 1980.g., u Šombierki rudniku uglja u basenu Gornja Silezija, predstavljala su rezultat normalnog rasedanja, sa različitim pružanjem ka severu i severozapadu (Slika 3-12).\r\nSlika 3-12. Važnije pojave rudarski generisanih potresa u Evropi i tipični fokalni mehanizmi (fokalni mehanizmi 1-4- [Gibowicz, 1984]; 5 - [Knoll, Kuhnt, 1990]; 6-7- [Hoang-Trang, i dr., 1988]; 8-9 - [Sileny, 1989] i 10 - [Hinzen, 1982]) [Wong, 1993; prilagođeno].\r\nPotres magnitude ML=4,6, koji se desio u Belčatov površinskom kopu 29. novembra 1980.g., nastao je kao rezultat reversnog ili moguće transkurentnog rasedanja, kao odgovor na opšte kompresione napone pravca severoistok-jugozapad (slika 12). Fokalni mehanizam drugog po veličini potresa u Šombierki rudniku, magnitude ML=4,1, koji je registrovan 12. jula 1981.g., pokazuje reversno rasedanje pružanja severoistok-\r\njugozapad, slično potresu u Belčatovu. Prema [Grunthal, Stromeyer, 1992], tektonsko naponsko polje u istočnoj Evropi karakteriše se kompresionim naponima u pravcu severoistok-jugozapad (krupne strelice na slici 3-12). Stoga, dva fokalna mehanizma (broj 1 i 4 na slici 3-12) ukazuju na to da su potresi u Belčatovu i Šombierki (1981) direktno pod uticajem tektonskog naponskog polja. S druge strane, normalna rasedanja za potrese broj 2 i 3 na slici 3-12, ukazuju na dominantnu ulogu litostatičkog pritiska. U svim slučajevima, međutim, s obzirom na to da se podzemne aktivnosti odvijaju na dubinama manjim od 1km, u radu [Gibowicz, 1984] predloženo je da tektonski naponi i u ovom slučaju moraju imati određenu ulogu, zbog male dubine, i, samim tim, malog litostatičkog pritiska. Sva četiri razmatrana potresa javila su se na dubini od 2km, ili manjoj, ukazujući na kompleksan odnos litostatičkog pritiska i tektonskih napona.\r\nOstale značajnije registrovane pojave rudarski generisanih potresa u Evropi uključuju (slika 3-12): (1) rudnike potaše u seizmički aktivnom južno-nemačkom bloku u Nemačkoj; (2) rudnike uglja u Rurskoj oblasti u zapadnoj Nemačkoj; (3) rudnike uglja u Provansi i basenu Lorene u Francuskoj; (4) rudnike uglja u oblasti Kladno, blizu Praga, i u oblasti Ostrava u severnoj Moravskoj (basen Gornja Silezija), u Češkoj. Svi registrovani potresi javljaju se u kompresionom naponskom polju orijentacije severoistok-jugozapad, u istočnoj Evropi, odnosno severozapad-jugoistok u zapadnoj Evropi [Miler, i dr., 1992].\r\nU južnonemačkom bloku, fokalni mehanizam razornih rudarskih potresa pokazuje reversno rasedanje sa pružanjem severozapad-jugoistok [Knoll, Kuhnt, 1990], što je u skladu sa orijentacijom regionalnog maksimalnog glavnog napona (Slika 12). Kompozitni fokalni mehanizam 21 malog potresa (ML<2,3) u Rurskoj oblasti takođe otkriva reversno rasedanje, ali kao odgovor na maksimalni glavni napon u pravcu severoistok-jugozapad [Hinzen, 1982]. Reversno rasedanje je takođe dobijeno fokalnim mehanizmima za potrese registrovane u rudnicima uglja u Provansi [Revalor, i dr., 1990], koje takođe odražava uticaj tektonskog naponskog polja. In situ merenja u ugljonosnom basenu u Provansi ukazuju na to da je maksimalni glavni napon oko 1,5 puta veći od litostatičkog pritiska [Gaviglio, i dr., 1990]. Najveći deo rudarski\r\ngenerisanih potresa u rudnicima uglja u basenu Lorene takođe je rezultat reversnog rasedanja [Hoang-Trang, i dr., 1988].\r\nAnaliza potresa u rudnicima Kladno i Ostrava Karvina u Češkoj [Sileny, 1989] ukazuje na predominantno normalno i transkurentno rasedanje. Seizmičnost registrovana u rudniku Kladno ne odražava uticaj tektonskog naponskog polja, s obzirom na to da normalno rasedanje ukazuje na ekstenzione uslove tektonskog polja. Relativno plitka žarišta ovih potresa (ispod 1km) mogu biti posledica relativnog malog uticaja tektonskih napona. S druge strane, potresi u Ostrava Karvina rudnicima pokazuju horizontalni položaj maksimalnog glavnog napona, što može odgovarati povećanom uticaju tektonskog polja na pojavu rudarski generisanih potresa (Slika 3- 12).\r\nŠto se tiče potresa u Južnoj Africi, kao oblasti sa najdubljom eksploatacijom mineralnih sirovina u svetu, u radu [Spottiswoode, i dr., 1971] autori su uočili da najveći deo registrovanih rudarski generisanih potresa u Vitvotersred oblasti predstavljaju rezultat loma na smicanje stenske mase. Ovi autori odredili su fokalne mehanizme za 11 potresa, magnitude ML=1-3,5, u Ist Rend rudnicima, i utvrdili normalna rasedanja, što se poklapa sa fokalnim mehanizmima tektonskih zemljotresa. U radu [Potgeiter, Roering, 1984] određeni su fokalni mehanizmi za 10 potresa, magnitude Ml=1,8-3,8, koji su se desili u periodu 1980-1981.g. u oblasti Klerksdorp. Fokalni mehanizmi pet potresa registrovani su u blizini postojećeg raseda i na udaljenosti 25m od iskopa, i pokazivali su normalno rasedanje, sa vertikalnim maksimalnim glavnim naponom. U radu [Brummer, Rorke, 1990] analizirano je pet velikih potresa, magnitude Ml=3,6-5,2, koji su registrovani u periodu 1982-1988.g. na području Južne Afrike. Svi potresi predstavljali su rezultat smicanja duž postojećih raseda, dok su četiri od pet potresa za koje je pomeranje duž raseda bilo moguće osmatrati, pokazivali normalno rasedanje. Prema tome, na osnovu dostupnih podataka, može se zaključiti da tektonsko naponsko polje nema dominantan uticaj na pojavu rudarskih potresa. To je pre svega posledica male seizmičke aktivnosti na području Južne Afrike, ali i velikih dubina eksploatacije, gde litostatistički pritisak ima dominantnu ulogu u formiranju naponskog polja.\r\n4. MODELOVANJE MEHANIZMA NASTANKA RUDARSKI GENERISANIH POTRESA DUŽ RASEDA\r\nDanas je opšte prihvaćeno stanovište da potresi duž raseda u Zemljinoj kori, tektonski ili indukovani, nastaju usled iznenadnog pada napona smicanja, koje je praćeno nestabilnim klizanjem duž raseda. U radu [Brace, Byerlee, 1966] prvi put je ukazano na ovu činjenicu, na osnovu rezultata laboratorijskog ispitivanja uzoraka granita u triaksijalnom aparatu. Mehanizam nastanka rupture, pri padu smičućeg napona, može se modelovati na dva načina: kao krti lom ili kao posledica trzanja usled naizmeničnog povećanja i smanjenja jačine trenja duž kontakta dva bloka stenske mase. Oba ova modela su matematički ekvivalentna, s obzirom na to da relativno kretanje povezuju sa padom napona smicanja duž rasedne površi. Međutim, sam proces nastanka pukotine razmatra se drugačije [Scholz, 2002]. U slučaju krtog loma, pretpostavlja se da je za nastanak pukotine neophodna određena karakteristična energija po jedinici površine, koja predstavlja materijalno svojstvo posmatranog sistema (stenske mase). U ovom slučaju, smatra se da stenska masa izvan ,,područja'' pukotine ostaje idealno elastična. Za slučaj idealno ravne, savršeno oštre pukotine nulte debljine, postoje tri osnovna načina kretanja duž kontakta dva dela stenske mase (slika 4-1). Tip 1 predstavlja tenzioni (otvarajući) tip, kod kojeg su pomeranja zidova pukotine upravna na pukotinu (slika 4-1a). Tip 2 predstavlja tip smicanja duž postojeće ravni, gde se pomeranja javljaju u ravni pukotine i normalna su na krajeve pukotina (slika 4-1b). Tip 3 predstavlja sličan slučaj sa prethodnim tipom, s tim što su sada pomeranja u ravni pukotine, paralelno njenim krajevima (slika 1c).\r\nS druge strane, u modelu trzajućeg kretanja1, pukotina se javlja kada napon duž rasedne površi dostigne vrednost koeficijenta statičkog trenja, pri čemu dolazi do pojave dinamičke nestabilnosti. U ovoj disertaciji, analitičko i numeričko modelovanje nastanka rasednih potresa vršiće se upravo pomoću modela koji opisuje ovo trzajuće kretanje.\r\nSlika 4-1. Tri osnovna načina nastanka pukotina u modelu krtog loma: (a) tenziono otvaranje; (b) smicanje u ravni pukotine upravno na njene krajeve; (c) smicanje u ravni pukotine, paralelno njenim krajevima [Gibowicz, Kijko, 1994].\r\n4.1. MODEL TRZANJA KAO MEHANIZAM NASTANKA RASEDNIH POTRESA\r\nU opštem slučaju, trzanje predstavlja kretanje koje se javlja na kontaktu dva čvrsta tela, usled naizmeničnog smanjenja i povećanja jačine trenja. Osnovni fenomenološki model kojim se opisuje trzajuće kretanje predstavlja Baridž-Knopof model bloka sa oprugom2 [Burridge, Knopoff, 1967]. U originalnom radu Baridža i Knopofa iz 1967.g., model se sastoji od osam blokova koji su međusobno povezani oprugama određene krutosti, i kreću se po hrapavoj površi, pri čemu je prvi blok, preko kotura povezan sa motorom, koji pokreće ceo sistem (slika 4-2).\r\nSlika 4-2. Shematski prikaz laboratorijskog BK modela; m predstavlja masu bloka, Fk - silu u opruzi, l - dužinu opruge, a V - brzinu kretanja sistema [Burridge, Knopoff, 1967].\r\nModel je konstruisan u cilju ispitivanja raspodele potencijalne energije u posmatranom sistemu, i poređenja te distribucije sa Gutenberg-Rihter3 i Omori-Utsu zakonom4 koji važe za tektonske zemljotrese. U tom smislu, razmatrana su dva slučaja: a) kada su blokovi povezani oprugama istih karakteristika; b) kada su blokovi povezani oprugama različitih karakteristika, pri čemu je uzeto da je krutost opruge proporcionalna zbiru masa svih blokova između opruge i slobodnog kraja. U datom opitu svi blokovi su imali masu od 142gr. Sila u oprugama u oba slučaja iznosila je 2N, pri čemu su u prvom opitu dužine svih opruga bile 3cm, dok su u drugom varirale od 1,5-12cm. Čitav sistem je povezan sa motorom koji uzrokuje njegovo kretanje. Brzina kretanja sistema iznosila je oko 2cm/min. U istom radu, Baridž i Knopof su odmah predložili i modifikaciju modela, tako da su umesto motora, sve blokove povezali oprugama određene krutosti za gornju pokretnu ploču (slika 4-3). Na taj način, postigli su kvalitativnu simulaciju kretanja na kontaktu tektonskih ploča ili duž rasedne zone. Upravo ovaj model se, u različitim varijantama, i danas koristi kao osnova u razmatranjima dinamike tektonskih zemljotresa.\r\nSlika 4-3. BK model sa pokretnom pločom (umesto motora) koja uzrokuje kretanje sistema, simulirajući kontakt dve tektonske ploče ili dva bloka stenske mase duž raseda.\r\nUkoliko usvojimo da je xim položaj i-tog bloka nakon kretanja čitavog sistema (tzv.,, lavine''), a h dužina opruge u nezategnutom stanju, onda je potencijalna energija sistema odmah nakon m-te lavine data sledećim izrazom:\r\nkoji ne predstavlja ništa drugo do primenu Hukovog zakona elastičnosti. Ovde x0 predstavlja položaj prvog bloka, koji je jedini izložen dejstvu spoljne sile motora Fp, i, kao takav, kreće se konstatnom brzinom tokom celokupnog trajanja eksperimenta. U tom smislu, može se smatrati da je proporcionalan vremenu.\r\nTokom trajanja eksperimenta, prvi blok deluje određenom silom F] na susedni blok preko opruge. Kada sila F] prevaziđe vrednost statičkog koeficijenta trenja, drugi blok počinje da se kreće, uzrokujući kretanje trećeg bloka, itd. Svaka promena energije nakon pomeranja blokova računa se po formuli (4.1). Shodno tome, moguće je odrediti količinu oslobođene energije u pojedinim vremenskim intervalima tokom kretanja sistema, kako je i prikazano na slici 4-4. U oba slučaja (kada su opruge jednake i različite dužine) jasno se uočava početni stadijum akumulacije potencijalne energije do trenutka kada pojedini blokovi počinju da se kreću. Međutim, za slučaj kada su opruge jednake dužine, vidimo da je početni stadijum akumulacije energije praćen kraćim intervalima oslobađanja energije (kada blokovi počinju da se kreću nezavisno) da bi od x0 = 120cm došlo do znatnog većeg oslobađanja akumulirane energije, što označava početak kretanja čitavog sistema. S druge strane, ukoliko su opruge nejednake dužine, jasno se uočava (slika 4-4b) da već od x0=40cm dolazi do kretanja celog sistema. Prema tome, može se zaklju čiti da je drugi slučaj znatno nepovoljniji po stabilnost sistema, i kao takav je često razmatran u savremenoj literaturi, najčešće putem variranja krutosti opruga.\r\nSlika 4-4. Potencijalna energija u funkciji pomeranja prvog bloka x0 za model više blokova povezanih oprugama, kada su: (a) sve opruge jednake dužine; (b) opruge različite dužine. [Burridge, Knopoff, 1967].\r\nTreba istaći da je pojava većih oslobađanja energije gotovo periodična, pri čemu se ta, gotovo, periodičnost može objasniti na sledeći način. Očigledno je da mora postojati određena gornja granica potencijalne energije, iznad koje je sistem uvek nestabilan, što je slučaj kada su sve opruge zategnute tako da su svi blokovi na granici stabilnosti i nestabilnosti. Stoga će bilo koji ,,okidač'' dovesti do naglog oslobađanja velike količine energije (tj. potresa u realnim uslovima). U prilog ovoj činjenici ide i grafik potencijalne energije sistema tokom dužeg vremenskog intervala, prikazan na slici 4-5, koji sa jedne strane podseća na ponašanje relaksacionog oscilatora, a s druge strane, na cikluse akumulacije i oslobađanja energije u blokovima stenske mase za vreme potresa.\r\nSlika 4-5. Skica oscilacija potencijalne energije u BK modelu u dužem vremenskom periodu.\r\nSada se logično nameće pitanje: na koji način su Baridž i Knopof povezali svoj model i mehanizam nastanka zemljotresa? Odgovor je prilično intiuitivan. Pretpostavljeno je da se količina oslobođene energije za vreme opita i magnituda potresa mogu smatrati analognim veličinama, s obzirom na to da postoji veliki broj empirijskih korelacija između ove dve veličine [Glavatović, 2001; Sunarić i Nedeljković, 2008], a potom je vršeno poređenje sa Gutenberg-Rihter zakonom, čime su dobijeni zadovoljavajući rezultati (Slika 4-6).\r\nSlika 4-6. Dijagram broja trzaja prema količini oslobođene energije za BK model, kada su: (a) sve opruge jednake dužine (R2=0,99); (b) opruge različite dužine (R2=0,98). Puna linija označava eksperimentalno dobijene podatke; isprekidanom linijom je predstavljena aproksimacija rezultata Gutenberg-Rihter zakonom. Razlika između krivih je uvećana zbog logaritamske razmere [Burridge, Knopoff, 1967; dopunjeno].\r\nNaime, u prvom slučaju, kada su opruge jednake dužine (slika 4-6a), dobijena je veoma dobra aproksimacija sa Gutenberg-Rihter zakonom (R2=0,99), sa vrednostima parametara, a=2,76 i 6=0,42, što je u skladu sa vrednostima za realno osmatrano tektonske zemljotrese. U drugom slučaju, kada su opruge nejednake dužine (slika 4-6b), kriva oslobođanja energije takođe je dobro aproksimirana Gutenberg-Rihter zakonom\r\n(R2=0,98), sa vrednostima parametara a=2,63 i &=0,18, koje su, takođe, u skladu sa vrednostima za realne zemljotrese.\r\nOsim poređenja sa Gutenberg-Rihter zakonom, Baridž i Knopof su takođe određivali količinu oslobođene potencijalne energije tokom usporavanja blokova do konačnog zaustavljanja sistema (slika 4-7a) za BK model sa 10 blokova i oprugama jednake dužine (iste krutosti). Dobijene rezultate poredili su sa Omori-Utsu zakonom opadanja energije naknadnih udara, čime su dobili zadovoljavajuću aproksimaciju, sa koeficijentom korelacije R2=0,92 (slika 4-7b).\r\nSlika 4-7. (a) Potencijalna energija u funkciji vremena određena tokom usporavanja blokova i konačnog zaustavljanja sistema. (b) Omori-Utsu zakon za ispitivani BK model. Na slici (b) puna linija označava eksperimentalno dobijene podatke; isprekidanom linijom je predstavljena aproksimacija rezultata Omori-Utsu zakonom [Burridge, Knopoff, 1967; dopunjeno].\r\nOve teorijsko-kvalitativne korelacije sa realno osmatranim zemljotresima, predložene u originalnom radu Baridža i Knopofa [1967], privukle su pažnju mnogih fizičara i matematičara i navele ih da studioznije pristupe ovom problemu. U daljem tekstu biće prikazani osnovi fenomenološko-matematičko-fizički modeli koji se danas koriste za izučavanje mehanizma nastanka zemljotresa.\r\nOvakvim mehanizmom reaktiviranja kretanja duž raseda mogu se opisati i rudarski generisani potresi. Pre početka iskopa, teren se nalazi u stanju ravnoteže. Nakon iskopa, stenska masa u zoni rasedanja počinje da se deformiše usled promene napona. Kako proces iskopa dalje napreduje, smanjuje se rastojanje između čela iskopa i raseda, uzrokujući dalju promenu napona. U određenom trenutku, stanje deformacije prelazi u stanje nestabilne ravnoteže, što dovodi do kretanja duž raseda i do pojave potresa. Tada blokovi stenske mase počinju da se kreću duž raseda trzanjem5, što je već predloženo kao način kretanja [Daihua, i dr. 1989; Lorig, Hobbs, 1990; Hsiung, i dr., 1990; McGarr, i dr. 2009]. Dalje kretanje duž raseda, pojava potresa, i njegovi efekti u velikoj i maloj razmeri vrlo su slični tektonskim zemljotresima [Kijko i dr., 1983]. Naime, dosadašnja istraživanja ukazuju na to da nema suštinskih razlika u mehanizmu nastanka između rudarski izazvanih potresa duž raseda i tektonskih zemljotresa, tako da se većina tehnika i postupaka iz seizmologije mogu primeniti pri proučavanju takvih, rudarski izazvanih potresa [Gibowicz, 1984; McGarr, 1984].\r\n4.2. TRENJE U STENSKOJ MASI KAO UZROK POJAVE NESTABILNOSTI\r\nBaridž i Knopof [1967] su pretpostavili da je čitava potencijalna nestabilnost sistema uslovljena trenjem između blokova i hrapave podloge. Pri tome, ovi autori su smatrali da je trenje funkcija samo brzine bloka u odnosu na hrapavu površ po kojoj se kreće (slika 4-8). Ova pretpostavka (da je trenje glavni izvor nelinearnosti i nestabilnosti sistema) predstavljaće osnovu čitavog jednog pravca daljeg istraživanja, čije ćemo rezultate delimično predstaviti u okviru ovog poglavlja.\r\nSlika 4-8. Skica promene trenja u funkciji brzine u BK modelu.\r\nZakon trenja zavisan od brzine, predložen u radu [Carlson, Langer, 1989] dat je u sledećem obliku:\r\nKao što se može videti, sila trenja se karakteriše sa dva parametra: S i a. Parametar S, uveden u radu [Carlson, 1991], predstavlja trenutni pad sile trenja na početku kretanja, dok parametar a opisuje stepen opadanja sile trenja pri porastu brzine kretanja6.\r\n4.3. ZAKONI TRENJA ZAVISNI OD BRZINE I STANJA\r\nDosadašnja razmatranja dinamike BK modela polazila su od činjenice da je trenje glavni izvor nelinearnosti sistema. Međutim, izraz za trenje dat jednačinom (4.5), ili prikazan na slici 4-8, i na prvi pogled deluje veoma jednostavno. S tim u vezi, nakon predloga modela bloka sa oprugom od strane Baridža i Knopofa, 1967.g., razvio se poseban pravac izučavanja prirode trenja između bloka i opruge. Ova istraživanja su naročito bila podstaknuta činjenicom da je bilo moguće izvoditi laboratorijske opite na uzorcima stenskih masa, i na taj način preciznije definisati zakone trenja. Upravo su eksperimentalna ispitivanja na uzorcima stenskih masa pokazala da trenje između bloka\r\ni hrapave podloge ne zavisi samo od brzine kretanja bloka, već i od stanja kontaktne površi. Rezultati ispitivanja prikazani su na slici 4-9.\r\nSlika 4-9. a) Promene vrednosti statičkog trenja u funkciji vremena kontakta dva bloka za inicijalno uglačane površine stenske mase (tamni simboli) i zdrobljeni kataklastični materijal (svetli simboli); b) Dijagram trenje-pomeranje koji ilustruje promenu koeficijenta statičkog trenja pri trzanju. Vreme statičkog kontakta (u stanju mirovanja blokova) izraženo je u sekundama, a brzina smicanja je 3p,m/s [Marone, 1998]\r\nNa slici 4-9a jasno se uočava zavisnost veličine koeficijenta statičkog trenja od vremena kontakta bloka i podloge. U sva četiri eksperimenta, za uzorke različitih vrsta stenskih masa, evidentno je da trenje raste tokom vremena, iako je blok u stanju mirovanja. Ova činjenica se objašnjava boljim prianjanjem kontakata između hrapavih površi7 blokova. Na slici 4-9b data je promena koeficijenta trenja tokom trajanja opita, kada se blokovi pomeraju trzajući. Jasno je da pri delovanju sile dolazi do naglog povećanja veličine trenja, i u momentu kada primenjena sila prevaziđe vrednost trenja, dolazi do kretanja blokova i naglog pada vrednosti trenja, koje sa usporavanjem blokova ponovo uspostavlja svoju prethodnu vrednost.\r\nPrethodno pomenuta osmatranja opisana su empirijskim konstitutivnim zakonom od strane Diteriha (1979) i Ruine (1983) u sledećem obliku:\r\nParametri A i B su empirijske konstante. Prema [Rice, 1993, 2001], parametar A predstavlja meru direktne zavisnosti od brzine (koja se ponekad naziva ,,direktni efekat''), dok je (A - B) mera zavisnosti trenja od brzine za vreme kretanja bloka. Parametar Dc predstavlja ,,kritično'' rastojanje - ono rastojanje koje je potrebno da blok pređe da bi se opet uspostavio statički kontakt hrapavih površi bloka i podloge. Parametar d predstavlja ,,promenljivu stanja'', koja oslikava uticaj hrapavosti podloge na kretanje bloka, i menja se na sledeći način:\r\nJednačina (4.7) zajedno sa jednačinom (4.6) predstavlja tzv. Diterihov zakon ili zakon usporavanja, zato što trenje raste i tokom stanja mirovanja, kada je V=0. Neophodno je međutim istaći da jednačina (4.6) nije definisana za slučaj kada je V=0. Shodno tome, da bi veličina trenja bila utvrđena, mora doći do bar malog kretanja duž posmatrane površi [Scholz i dr., 1972; Baumberger i dr., 1994]. Na slici 4-10 data je grafička ilustracija Diterihovog zakona trenja.\r\nU radu [Ruina, 1983] predložen je malo drugačiji zakon trenja, za koji su brzina i pomeranje od primarnog značaja (za razliku od Diterihovog zakona gde je vreme imalo primarnu ulogu):\r\nU ovom slučaju, Ruina smatra da bilo koja promena trenja zahteva kretanje duž površi. Drugim rečima, trenje se ne menja za vreme stanja mirovanja bloka, kao u Diterihovom zakonu trenja.\r\nSlika 4-10. Dijagram zavisnosti koeficijenta trenja ^ od pomeranja u i brzine kao ilustracija odgovora sistema na iznenadno povećanje brzine smicanja. Primenjena brzina, inicijalno održavana konstantnom, V0, iznenada se povećava za AV, i potom se održava konstantnom V0 + AV. Koeficijent trenja, inicijalno konstantan, T0, iznenada se poveća za vrednost A, a potom eksponencijalno opada, uzimajući novu vrednost B. Dužina Dc predstavlja rastojanje potrebno da se vrednost promenljive stanja promeni od do do 0.\r\nMeđutim, iako je razlika između zakona trenja Diteriha [1979] i Ruine [1983] fundamentalna u smislu mikromehaničke interpretacije procesa, oba zakona ,,reprodukuju'' laboratorijske podatke na sličan način (slika 4-11). Drugim rečima, u oba slučaja zapaža se tzv. direktni efekat, kada, usled povećanja sile koja deluje na blok, raste i koeficijent trenja do neke vršne vrednosti, nakon čega dolazi do naglog pada vrednosti trenja, koje dostiže minimum kada se blok konačno zaustavi, a potom se ceo proces ponavlja. Zakoni se jedino razlikuju u karakteru i veličini promena koeficijenta trenja. Naime, kod Diterihovog zakona trenja, imajući u vidu značaj vremena kontakta između bloka i podloge, kriva promene trenja u zavisnosti od pomeranja je asimetrična (prva kriva na slici 4-11). S druge strane, kod Ruina zakona trenja, povećanje vrednosti statičkog trenja u stacionarnom stanju je nezavisno od vremena, i, stoga je simetrično u odnosu na promene brzine (druga kriva na slici 4-11).\r\nSlika 4-11. Trenje u zavisnosti od normalizovanog pomeranja za tri zakona trenja zavisna od brzine i stanja. Konstitutivni parametri definisani su na slici za prvu krivu, i važe i za ostale razmatrane slučajeve. Za početnu brzinu kretanja usvojeno je da je Vo=1pm/s. Krive su dobijene za sledeće vrednosti parametara: A=0,01; 5=0,015; Dc=20^m, i *=0,01p/m.\r\nTreći zakon trenja predložen je u radu [Perrin, Rice, Zheng, 1995], u sledećem obliku:\r\nOvaj zakon trenja daje slične rezultate kao i drugi zakoni (slika 4-11), sa opadanjem vrednosti trenja do stacionarnog stanja po zakonu e(~u/Dc), gde je u pomeranje.\r\nU ovoj disertaciji, pri razmatranju BK modela sa zakonom trenja zavisnim od brzine i stanja, usvaja se Diterihov zakon trenja, koji se vrlo često naziva i Diterih-Ruina zakonom trenja.\r\nNeophodno je napomenuti da Diterihov zakon trenja opisuje promenu koeficijenta trenja ne samo za površi stenskih masa, već i za metale [Popov i dr., 2010], listove papira [Heslot, i dr., 1994], i dr.\r\nMeđutim, i pored široke univezalnosti i aplikativnosti ovog zakona ponašanja materijala8, jedan broj eksperimenata ukazuje na činjenicu da Diterihov zakon trenja ne važi u uslovima velike brzine klizanja, što se objašnjava različitim mehaničko- hemijskim reakcijama (zagrevanje tela9, frikciono topljenje10, stvaranje silikatne želatinske mase duž površi kontakta11, stvaranje granularnog sloja na kontaktu, koji odgovara tektonski zdrobljenim stenskim masama12) koje nastaju usled frikcionog zagrevanja, koje ,,podmazuje'' kontaktne površine do nekog stepena (Di Toro i dr., 2004; Mizoguchi, i dr.., 2006;). Pojedini autori smatraju da je uloga frikcionog zagrevanja od velikog značaja pri kretanju duž raseda, s obzirom na to da su pritisci u seizmogenoj zoni reda veličine 100MPa. Međutim, naše znanje o ovim pojavama još uvek nije potpuno, tako da se smatra da je, u laboratorijskim uslovima, uzimajući u obzir razmeru posmatranja, Diterih-Ruina zakonom trenja moguće fenomenološki opisati glavni mehanizam nastanka potresa.\r\nPored Baridž-Knopof modela bloka sa oprugom, vrlo često se koriste i tzv. mrežni modeli sa spregnutim preslikavanjima13, poznatiji kao modeli samo-organizovane kritičnosti14 (Bak, Tang, Wiesenfeld, 1987; Nakanishi, 1990; Olami, Feder, Christensen, 1992; Hergarten, Neugebauer, 2000; Helmstetter, Hergarten, Sornette, 2004), sa najpoznatijim OFC modelom, potom model svežnja vlakana15 (Pradhan, i dr., 2010;\r\nHemmer, Hansen, 1992), kao i model preklapanja dve fraktalne površi16 (Bhattacharya, Chakrabarti, 2006; Pradhan, i dr., 2003; Kawamura i dr., 2012). Međutim ovi modeli ne opisuju trzajuće kretanje pri nastanku zemljotresa, tako da ne mogu služiti za modelovanje rasednih potresa razmatranih u ovoj disertaciji.\r\nU ovoj disertaciji vrši se razmatranje Baridž-Knopof modela sa jednim i dva bloka, kod kojeg je trenje opisano sa dva zakona: (a) zakonom zavisnim od brzine, kao u originalnom radu Baridža i Knopofa, i (b) Diterih-Ruina zakonom trenja.\r\n5. NELINEARNA DINAMIKA I TEORIJA HAOSA - PREGLED OSNOVNIH POJMOVA\r\nS obzirom na to da su u ovoj disertaciji prvi put primenjene metode nelinearne dinamike u rudarskom inženjerstvu, u ovom poglavlju je dat kratak prikaz najvažnijih pojmova iz oblasti nelinearne dinamike i teorije haosa [Strogatz, 1994; Sprott, 2003].\r\nNELINEARNI DINAMIČKI SISTEM - nelinearni dinamički sistem predstavlja dinamički sistem čije se ponašanje opisuje skupom nelinearnih jednačina. Drugim rečima, dinamičke promenljive koje opisuju svojstva sistema (položaj, brzina, ubrzanje, pritisak, i dr.) opisuju se jednačinama u nelinearnom obliku. Pri tome, postoje dva osnovna tipa dinamičkih sistema:\r\nO diferencijalne jednačine, koje opisuju ponašanje sistema u kontinualnom vremenu;\r\nO diferencne jednačine (iterativna preslikavanja), koje opisuju ponašanje sistema u diskretnim vremenskim intervalima. Pored ove podele, u odnosu na karakter vremenskog intervala posmatranja sistema, i diferencijalne i diferencne jednačine mogu biti:\r\nO autonomne, kada ne pokazuju eksplicitnu vremensku zavisnost; O neautonomne, kod kojih ponašanje sistema eksplicitno zavisi od vremena. U ovoj disertaciji, mehanizam nastanka rudarski generisanih potresa opisan je sistemima običnih autonomnih diferencijalnih jednačina i diferencijalnih jednačina sa kašnjenjem prvog reda.\r\nBIFURKACIJE - bifurkacije predstavljaju kvalitativne, topološke promene ponašanja dinamičkih sistema, pri promeni vrednosti kontrolnih parametara (koje se nazivaju bifurkacionim vrednostima) za vrednost u samoj bifurkaciji. Broj parametara, čijom promenom dolazi do pojave bifurkacija, predstavlja kodimenziju bifurkacije. Među bifurkacijama kodimenzije 1, kao najčešće, izdvajaju se sedlo-čvor bifurkacija, transkritična i vilasta bifurkacija, Hopfova bifurkacija i bifurkacija sa udvajanjem perioda, dok među bifurkacijama kodimezije 2 obično razlikujemo homokliničku i heterokliničku bifurkaciju. Strogo matematički posmatrano, do pojave bifurkacija dolazi kada sopstvena vrednost karakteristične jednačine posmatranog sistema menja predznak (prolazi kroz nultu vrednost). U slučaju konjugovano-kompleksnih rešenja karakteristične jednačine, realni deo, Re(X) mora da promeni predznak da bi došlo do pojave bifurkacije.\r\nKarakteristični primeri lokalnih bifurkacija.\r\nO sedlo-čvor bifurkacija - kod ovog tipa bifurkacije dolazi do pojave i nestanka fiksnih tačaka, pri promeni vrednosti kontrolnog parametra (slika 5-1a); O transkritična bifurkacija - u ovom slučaju dolazi do pojave, spajanja i ponovnog razdvajanja fiksnih tačaka (slika 5-1b);\r\nSlika 5-1. (a) sedlo-čvor bifurkacija; (b) transkritična bifurkacija.\r\nO vilasta bifurkacija - u ovom slučaju razlikujemo:\r\nO natkritičnu vilastu bifurkaciju, kada sa promenom kontrolnog parametra fiksna tačka gubi stabilnost, uz formiranje dve nove stabilne fiksne tačke (slika 5-2a);\r\nO potkritična vilasta bifurkacija, kada fiksna tačka gubi stabilnost, uz nestajanje dve nestabilne fiksne tačke (slika 5-2b).\r\nSlika 5-2. (a) natkritična vilasta bifurkacija; (b) potkritična vilasta bifurkacija.\r\nU ovoj disertaciji razmatraju se lokalne bifurkacije (Hopfova bifurkacija, u prvom redu) koje se javljaju pri promeni kontrolnih parametara posmatranih matematičkih modela nastanka rudarski generisanih potresa.\r\nFIKSNA (SINGULARNA) TAČKA - u matematici, tačka x€Rn predstavlja fiksnu (ravnotežnu) tačku za diferencijalnu jednačinu dx/dt=f(x), ukoliko je f(x)=0 za svako t. U fiksnoj tački, sistem se nalazi u ravnotežnom stanju. Obično se fiksne tačke klasifikuju prema predznaku sopstvenih (karakterističnih) vrednosti linearizovanog sistema diferencijalnih jednačina u blizini ravnotežnog stanja. Drugim rečima, stabilnost fiksne tačke se određuje nalaženjem Jakobijan matrice sistema u fiksnoj tački, a potom određivanjem njenih karakterističnih (sopstvenih) vrednosti. Fisknu tačku nazivamo hiperboličkom ukoliko su realni delovi sopstvenih rešenja različiti od nule:\r\nO ukoliko realni delovi svih sopstvenih rešenja imaju negativan predznak, onda kažemo da je fiksna tačka nestabilna (izvor ili nestabilni fokus); O ukoliko realni delovi svih sopstvenih rešenja imaju pozitivan predznak, onda kažemo da je fiksna tačka stabilna (atraktor ili stabilni fokus). Klasifikacija fiksnih tačaka ima vrlo jednostavnu grafičku interpretaciju. Naime, posmatrajmo jedan linearni sistem zapisan u matričnom obliku:\r\nOdgovarajuća karakteristična jednačina sistema glasi:\r\ngde je T=a+d, a A=ad- bc. Sopstvene vrednosti (rešenja karakteristične jednačine) u ovom slučaju imaju sledeći oblik:\r\nLako se uočava da su u (T, A) ravni parovi realnih i konjugovano-kompleksnih sopstvenih vrednosti odvojeni parabolom T2=4A, a da se granični slučajevi mogu pojaviti kada je A=0, odnosno T=0 (slika 5-3).\r\nSlika 5-3. Grafički prikaz klasifikacije fiksnih tačaka.\r\nU zavisnosti od predznaka sopstvenih vrednosti, pojavljuju se različiti tipovi fiksnih tačaka:\r\n• Realne sopstvene vrednosti: T2>4A, kada se sopstvene vrednosti nalaze ispod parabole u (T,A) ravni. Tada se mogu razlikovati dva karakteristična slučaja:\r\n- Slučaj I: A>0, kada se može lako pokazati da su sopstvene vrednosti istog znaka, jer je:\r\n1. ako je T<0, fiksna tačka je stabilni čvor (slika 5-4a);\r\n2. ako je T>0, fiksna tačka je nestabilni čvor (slika 5-4b). Geometrijsko mesto ovih tačka ograničeno je parabolom T2=4A i pravom A=0.\r\nSlika 5-4. Tip fiksne tačke za A=0, kada su sopstvene vrednosti istog znaka: (a) stabilni čvor; (b) nestabilni čvor.\r\n- Slučaj II: A<0. U ovom slučaju sopstvene vrednosti su različitog znaka jer je:\r\nTada stacionarna tačka predstavlja tačku-sedlo bez obzira na znak parametra T (slika 5-5). Geometrijsko mesto ovih fiksnih tačaka ograničeno je pravom A=0, i obuhvata celu poluravan određenu sa A<0.\r\nSlika 5-5. Tip fiksne tačke za A<0, kada su sopstvene vrednosti različitog znaka (tačka- sedlo).\r\n- Granični slučajevi. Ako se analiziraju realne sopstvene vrednosti, mogu se uočiti sledeći granični slučajevi:\r\n■ Slučaj I: T2=4A. U ovom slučaju dobijaju se parovi dvostrukih realnih sopstvenih vrednosti l]=l2=r/2, koje odgovaraju degenerisanim čvorovima, koji su stabilni kada je T<0, a nestabilni kada je T>0 (slika 5-6). Geometrijsko mesto ovih fiksnih tačaka je parabola T2 = 4A.\r\nSlika 5-6. Tip fiksne tačke za T2=4A, za parove dvostrukih realnih sopstvenih vrednosti: (a) stabilni degenerisani čvor (T<0); (b) nestabilni degenerisani čvor (T>0).\r\n■ Slučaj II: A=0. U ovom slučaju, imamo sledeće sopstvene vrednosti: Aj=0 i X2=T, koje odgovaraju neizolovanim čvorovima koji su stabilni kada je T<0, a nestabilni kada je T>0 (slika 5-7).\r\nSlika 5-7. Tip fiksne tačke za A=0, kada je X=0 i X2=T\\ (a) stabilni neizolovani čvorovi (t<0); (b) nestabilni neizolovani čvorovi (t>0).\r\n■ Slučaj III: A<0, T=0. U ovom slučaju dobija se par sopstvenih vrednosti X12=± A , koje su realne i različite i odgovaraju tački- sedlu. Dakle, oblast ispod apscise u kojoj se nalazi ovaj tip stacionarnih tačaka ostaje nepodeljena. Ako se analiziraju konjugovano-kompleksne sopstvene vrednosti, granični slučaj se dobija za T=0. Tada imamo par imaginarnih sopstvenih vrednosti X1j2=± iA , koje odgovaraju ravnotežnim tačkama tipa centra (slika 5-8). Njihovo geometrijsko mesto je pozitivni deo ordinate A >0, T=0.\r\nSlika 5-8. Tip fiksne tačke za T=0, kada su sopstvene vrednosti XIT2=± iA (fiksna tačka tipa centra).\r\n• Konjugovano-kompleksne vrednosti: T2<4A, kada se sopstvene vrednosti nalaze iznad parabole u (T,A) ravni, čime je istovremeno zadovoljena i nejednakost A>0. U ovom slučaju, fiksne tačke se, s obzirom na znak parametra t, mogu podeliti u dve grupe:\r\n1. ako je T<0, fiksna tačka je stabilni fokus (slika 5-9a);\r\n2. ako je T>0, fiksna tačka je nestabilni fokus (slika 5-9b). Geometrijsko mesto ove dve grupe sopstvenih vrednosti ograničeno je parabolom T2=4A i pravom A=0.\r\nSlika 5-9. Tip fiksne tačke za T2<4A: (a) stabilni fokus (t<0); (b) nestabilni fokus (t>0).\r\nU ovoj disertaciji, fiksne tačke se klasifikuju samo kao stabilne ili nestabilne, dok detaljnija klasifikacija, kao na slici 5-3, neće biti vršena.\r\nGRANIČNI CIKLUS - granični ciklus predstavlja izolovanu zatvorenu trajektoriju, ka kojoj se susedne trajektorije ili približavaju (stabilni granični ciklus) ili se od nje udaljavaju (nestabilni granični ciklus). Granični ciklusi su karakteristična odlika nelinearnih sistema - ne mogu se javiti u linearnim. U ovoj disertaciji, granični ciklusi neće biti posebno razmatrani, s obzirom na to da predstavljaju promenu ponašanja sistema na većoj udaljenosti od fiksne tačke (dok je, u ovom slučaju, analiza dinamike ograničena samo na linearizaciju sistema i ispitivanju njegovog ponašanja za malu okolinu fiksne tačke, tzv. lokalnu analizu). Definicija graničnih ciklusa data je u ovom poglavlju, jer se javljaju pri prolasku sistema kroz Hopfovu bifurkaciju.\r\nFAZNI PROSTOR - fazni prostor predstavlja višedimenzionalni prostor u kome su predstavljena sva moguća stanja sistema.\r\nHOPFOVA BIFURKACIJA - Hopfova bifurkacija ili Andronov-Hopfova bifurakcija predstavlja lokalnu bifurkaciju kada fiksna tačka dinamičkog sistema gubi stabilnost i kada par konjugovano-kompleksnih sopstvenih vrednosti linearizovanog sistema u okolini fiksne tačke prelazi preko imaginarne ose kompleksne ravni.\r\nU opštem slučaju, pođimo od obične diferencijalne jednačine prvog reda:\r\nkoja zavisi od parametra aCR, gde je f glatka funkcija. Neka je x0(a) fiksna (ravnotežna) tačka sistema. Pretpostavimo da Jakobijan matrica J(a)=fx(x0(a),a) ima jedan par konjugovano-kompleksnih svojstvenih vrednosti:\r\nKada a menja predznak (prolazi kroz tačku a=0), ravnotežno stanje sistema menja stabilnost, i dolazi do pojave graničnog ciklusa.\r\nU opštem slučaju, postoje dva osnovna tipa Hopf bifurkacije - natkritična i potkritična. Posmatrajmo sistem običnih diferencijalnih jednačina prvog reda, koje zavise od parametra p.:\r\nPretpostavimo da je fiksna tačka sistema u koordinatnom početku za sve vrednosti Pretpostavimo takođe da Jakobijan matrica sistema Df ima dve imaginarne sopstvene vrednosti Ai(ju) i A2(JU) kada je p.=p.c. Ukoliko su realni delovi sopstvenih vrednosti pozitivni, i ukoliko je koordinatni početak asimptotski stabilna fiksna tačka, p.=p.c, onda je:\r\nO p.c bifurkaciona tačka;\r\nO za neko takvo da je < p. < p.c, koordinatni početak je stabilni fokus;\r\nO za neko takvo da je p.c < n < p2, koordinatni početak je nestabilan, okružen stabilnim graničnim ciklusom, čija se veličina povećava sa povećanjem vrednosti\r\nU ovom slučaju, Hopfova bifurkacija je natkritična (slika 5-10a).\r\nS druge strane, ukoliko su realni delovi svojstvenih vrednosti pozitivni, i ukoliko je koordinatni početak asimptotski stabilna fiksna tačka, p.=p.c, onda je:\r\nO p.c bifurkaciona tačka;\r\nO za neko takvo da je ^ < n < p.c, koordinatni početak je stabilni fokus, s tim da na većoj udaljenosti od koordinatnog početka postoje i nestabilni i stabilni granični ciklus;\r\nO za neko takvo da je p.c < n < p2, koordinatni početak je nestabilan, okružen stabilnim graničnim ciklusom, čija se veličina povećava sa povećanjem vrednosti /u, i koji je na znatno većem rastojanju od novostvorene nestabilne fiksne tačke, nego u slučaju natkritične Hopfove bifurkacije.\r\nU ovom slučaju, Hopfova bifurkacija je potkritična (slika 5-10b)\r\nSlika 5-10. Hopfova bifurkacija: (a) natkritična; (b) potkritična.\r\nTORUS - torus predstavlja površ koja nastaje neprekidnim rotacijama kruga k, poluprečnika r, oko jedne ose koja se nalazi u njegovoj ravni (slika 5-11). U parametarskom obliku, jednačine torusa su:\r\ngde su:\r\n- - uglovi koji prave pun krug, sa početkom u 0, i krajem u 2n, tako da njihove vrednosti počinju i završavaju se u istoj tački;\r\n- R - rastojanje od centra ,,cevi''do centra torusa (veliki poluprečnik torusa);\r\n- r - radijus cevi (mali poluprečnik torusa).\r\nSlika 5-11. Ilustracija nastanka torusa.\r\nDETERMINISTIČKI HAOS - iako ne postoji opšte prihvaćena definicija termina haos, većina istraživača je prihvatila sledeću grubu definiciju: ,,Deterministički haos predstavlja aperiodično dugoročno ponašanje u determinističkom sistemu, koje je osetljivo na male promene početnih uslova, i koje, u konačnom vremenskom intervalu konvergira ka stranom atraktoru''. Drugim rečima, za razliku od opšte prihvaćenog termina (razvijeni) haos, kao sinonima za potpunu neuređenost (maksimalnu entropiju) sistema, deterministički haos predstavlja ponašanje koje ima određeni mehanizam ,,u pozadini'' i koje je ograničeno na nekom konačnom skupu. Da bi došlo do pojave haotičnog ponašanja, sistem mora da prođe kroz niz bifurkacija. Lokalnim bifurkacijama se ,,do haosa'' dolazi najčešće na tri osnovna načina: bifurkacijama sa udvajanjem perioda, kvaziperiodičnim i intermitentnim putem [Sprott, 2003]. Bifurkacije sa udvajanjem perioda javljaju se kod logističkog preslikavanja, gde stabilne fiksne tačke postaju nestabilne nakon niza flip bifurkacija (bifurkacija kroz koju fiksna tačka gubi stabilnost i dolazi do nestanka graničnog ciklusa, sa promenom vrednosti parametara). Kvaziperiodičnost se javlja kod Van der Polovog oscilatora, gde spiralni čvor postaje nestabilan i javlja se granični krug kroz Hopfovu bifurkaciju [Marsden, McCracken, 1976], koji potom prelazi u torus kroz drugu Hopfovu bifurkaciju, nakon čega dolazi do raspadanja torusa i pojave deterministički haotičnog ponašanja [Ruelle, Takens, 1971; Newhouse, i dr., 1978]. Ovaj ,,put u haos'' je poznat kao Ruel-Takens-\r\nNjuhaus put1 u haos, i pojava determinističkog haosa u modelima potresa istraživanim u ovoj disertaciji upravo će se odigravati na ovaj način (slika 5-12).\r\nSlika 5-12. Ruel-Takens-Njuhaus put u haos.\r\nTreći, intermitentni put u haos, karakteriše se sedlo-čvor bifurkacijom u višedimenzionalnim sistemima, pri kojoj se sedlo tačka i stabilni čvor spajaju i nestaju, dajući trajektoriju sa ,,haotičnim'' periodima koja je u smeni sa regularnim oscilacijama [Pomeau, Manneville, 1980].\r\nU ovoj disertaciji, u slučaju izučavanja Madariaga sistema sa jednim blokom, put u haos je predstavljen upravo Ruel-Takens-Njuhaus „scenarijom\". Kod proučavanja drugog modela, sa jednim i dva bloka, kao i trećeg modela sa jednim blokom put u haos se posebno ne definiše, već se samo osmatra vrednost kontrolnih parametara za koje dolazi do pojave determinističkog haotičnog ponašanja.\r\nSTRANI ATRAKTOR - strani atraktor predstavlja skup konačne veličine ka kojem sve susedne trajektorije (rešenja) sistema konvergiraju. Strogo posmatrano, strani atraktor se definiše kao zatvoreni skup A, sa sledećim svojstvima:\r\nO A je invarijantni skup: bilo koja trajektorija x(t) koja počinje u A ostaje u A sve vreme;\r\nO A privlači otvoreni skup početnih uslova: postoji otvoreni skup U koji sadrži A tako da ako je x(0)CU, onda rastojanje od x(t) do A teži nuli kada t ^ To znači da A privlači sve trajektorije koje počinju dovoljno blizu skupu A. Najveći takav skup U naziva se bazom atrakcije skupa A; O Skup A je minimalan: ne postoji odgovarajući podskup skupa A koji zadovoljava uslove 1 i 2;\r\nO Skup A pokazuje osetljivost na malu promenu početnih uslova; O Skup A ima fraktalnu dimenziju. Strani atraktor takođe može biti rekonstruisan iz skalarne vremenske serije, o čemu će biti više reči u poglavlju ,,Metodologija istraživanja''. U ovoj disertaciji, neće biti određivana vrsta ni dimenzija stranog atraktora, već je primarni naglasak stavljen na pojavu determinističkog haosa.\r\nLJAPUNOVLJEV EKSPONENT - Ljapunovljev eksponent predstavlja veličinu eksponencijalne divergencije dve susedne trajektorije sistema. Pretpostavimo da je x(t) tačka na stranom atraktoru u vremenskom trenutku t. Označimo sa x(t)+6(t) početno blisku (susednu) tačku, gde 6 predstavlja vektor rastojanja, početne dužine reda veličine =1015 (slika 5-13).\r\nSlika 5-13. Eksponencijalna divergencija između dve susedne trajektorije (rešenja) sistema.\r\nUkoliko je atraktor strani (determinističko haotično ponašanje), promena rastojanja S(t) tokom vremena može se opisati na sledeći način: | D(t) |~|S0| eAt, gde X predstavlja Ljapunovljev eksponent sistema. Zapravo, za «-dimenzionalni sistem postoji n različitih vrednosti Ljapunovljevog eksponenta. Naime, ukoliko razmatramo ponašanje infinitezimalno male sfere početnih uslova, tokom vremena će doći do deformacije početne sfere u infinitezimalni elipsoid. Neka Sk(t), k=1,...,n označava dužinu k-te ose elipsoida. Onda je Sk(t)~Sk(0)eAkt, gde su Xk Ljapunovljevi eksponenti. Za veliko t, prečnik elipsoida je određen najvećom pozitivnom vrednošću Xk. Stoga, traženo X zapravo predstavlja najveći (maksimalni) Ljapunovljev eksponent. Neophodno je takođe napomenuti da kada sistem ima pozitivan Ljapunovljev eksponent, postoji određeni vremenski interval nakon koga predviđanje ponašanja sistema više nije moguće. Ukoliko prepostavimo početno malu grešku u merenju determinističkog haotičnog sistema, tokom vremena ta greška će eksponencijalno brzo da raste (| $(t) I ~ I So I eAt). Ukoliko je a mera naše tolerancije, tada smatramo da kada je | S(t) | > a, onda predviđanje ponašanja sistema postaje nemoguće, što se događa nakon vremenskog intervala reda veličine:\r\nVrednosti maksimalnog Ljapunovljevog eksponenta kvalitativno ukazuju na ponašanje posmatranog sistema:\r\n■ X < 0 - trajektorije (rešenja) sistema teže stabilnoj fiksnoj tački ili stabilnom graničnom ciklusu (slika 5-14). Negativni Ljapunovljevi eksponenti su karakteristični za disipativne ili ne-konzervativne sisteme (prigušeni harmonijski oscilator). Takvi sistemi pokazuju asimptotsku stabilnost - što je ,,negativnija'' vrednost eksponenta, veća je stabilnost sistema. Superstabilne fiksne tačke i superstabilni granični ciklusi imaju Ljapunovljev eksponent - ro. Ova vrednost se, na primer, javlja kod prigušenog harmonijskog oscilatora koji teži ravnotežnom stanju najvećom mogućom brzinom.\r\nSlika 5-14. Skica trajektorija (rešenja) sistema za koje je X < 0; a) stabilna fiksna tačka (spirala); b) stabilni granični ciklus.\r\n■ X = 0 - trajektorija predstavlja neutralnu fiksnu tačku - centar (slika 5-15). Nulta vrednost Ljapunovljevog eksponenta ukazuje na to da je sistem u nekoj vrsti stacionarnog stanja. Fizički sistem sa nultom vrednošću ovog eksponenta je konzervativan. Ukoliko uzmemo primer dva identična harmonijska oscilatora sa različitim amplitudama, onda će njihovi fazni portreti biti predstavljeni u vidu koncentričnih krugova.\r\nSlika 5-15. Skica trajektorija rešenja sistema za koje je X = 0 - neutralna fiksna tačka (centar).\r\n■ X > 0 - trajektorija rešenja sistema je nestabilna. Susedne vrednosti (tačke) ma koliko bile bliske, tokom vremena će pokazivati divergenciju, odnosno determinističko haotično ponašanje (slika 5-16).\r\nSlika 5-16. Skica trajektorije rešenja za deterministički haotičan sistem (Lorencov sistem), za koji je X > 0 (strani atraktor).\r\nU ovoj disertaciji, Ljapunovljev eksponent se određuje pomoću dve metode - Vulfove i Rozenštajnove, koje su detaljno opisane u poglavlju ,,Metodologija istraživanja''.\r\nLORENCOV SISTEM - Lorencov sistem je sistem običnih diferencijalnih jednačina, predložen od strane Edvarda Lorenca, 1963.g., koji je, proučavajući pojednostavljeni model konvekcionih strujanja u atmosferi, došao do zaključka da rešenja njegovih jednačina nikada ne konvergiraju ravnotežnom ili periodičnom stanju - umesto toga, rešenja pokazuju dugotrajne, neregularne, aperiodične oscilacije. Kada strani atraktor Lorencovog sistema projektujemo na dvodimenzionalnu ravan, trajektorije (rešenja) podsećaju na krila leptira2. Danas se Lorencov sistem smatra ,,svetim gralom'' nelinearne dinamike i teorije haosa, a Edvard Lorenc jednim od najvažnijih istraživača- pionira u ovoj oblasti.\r\nLorencov sistem je dat u sledećem obliku:\r\ngde su a, r i fi parametri sistema, pri čemu je Lorenc utvrdio da za vrednosti parametara o=10, r=28 i fi=8/3 ispitivani sistem pokazuje deterministički haotično ponašanje. Strani atraktor sistema prikazan je u faznom prostoru na slici 5-17.\r\nSlika 5-17. Strani atraktor Lorencovog sistema (,,krila leptira'').\r\nDIFERENCIJALNE JEDNAČINE SA KASNJENJEM3 - diferencijalne jednačine sa kašnjenjem pripadaju grupi funkcionalnih diferencijalnih jednačina, opšteg oblika:\r\ngde je x(t) € Rn, i f predstavlja nelinearnu glatku funkciju, koja zavisi od (m+1) diskretnih vrednosti promenljivih x, kao i od autonomnih (vremenski nezavisnih) parametara v € R. Vremenska kašnjenja r(x(t)) mogu biti ili konstante ili kašnjenja zavisna od stanja. U ovoj disertaciji, proučavane diferencijalne jednačine imaju konstanto kašnjenje T.\r\nZa razliku od običnih diferencijalnih jednačina, gde je jedinstveno rešenje određeno vrednošću jedne tačke u Euklidskom prostoru za početni trenutak t0, za rešavanje diferencijalnih jednačina sa kašnjenjem, kao početni uslov potrebni su podaci o celom intervalu [t0-T,t0]. Drugim rečima, da bismo znali vrednost izvoda u trenutku t0, mora nam biti poznata vrednost početne funkcije na celom intervalu [t0-T, t0]. Shodno tome, da bi početni uslov u slučaju rešavanja diferencijalnih jednačina sa kašnjenjem imao smisla, potrebno je da početna funkcija ili početna istorija budu poznati - vrednost x(t0) u intervalu [t0-T,t0]. Svaka takva početna funkcija određuje jedinstveno rešenje diferencijalne jednačine sa kašnjenjem.\r\n6. METODOLOGIJA ISTRAZIVANJA\r\nMetode korišćene u ovoj disertaciji mogu se svrstati u dve klase. U prvu klasu metoda spada matematičko modelovanje fenomenološkog modela mehanizma nastanka potresa (Baridž-Knopof modela), koje se operativno sprovodi analitičkim i numeričkim rešavanjem sistema običnih diferencijalnih jednačina, sa i bez kašnjenja, koje opisuju kretanje bloka. Drugim rečima, vrši se analiza ponašanja rešenja jednačina blizu ravnotežnog stanja (numerički), kao i standardna lokalna bifurkaciona analiza (analitički), dalje verifikovana numeričkim putem korišćenjem DDE-BIFTOOL-a. Pri tome, kvalitativne promene dinamike modela (bifurkacije) razmatraju se ili usled periodičnih perturbacija odabranih kontrolnih parametara, ili kao posledica uvođenja novih parametara (jačina trenja c i vremensko kašnjenje T).\r\nDruga klasa metoda pripada statističkim metodama. Primenom ovih metoda vrši se izučavanje registrovanih seizmograma oscilovanja tla za vreme rudarski generisanih potresa. Konkretno, primenjuju se dva različita postupka - analiza nelinearnih vremenskih serija i analiza surogat podataka.\r\nU prilogu br.2 dat je algoritam metodologije nelinearnog dinamičkog modelovanja potresa izazvanih promenom naponskog stanja pri izradi horizontalnih podzemnih prostorija.\r\nU ovom poglavlju detaljno se opisane opšte metode ispitivanja, imajući u vidu činjenicu da se metode nelinearne dinamike prvi put praktično primenjuju u rudarskom inženjerstvu.\r\n6.1. ANALIZA REGISTROVANIH OCILOVANJA TLA ZA VREME POTRESA 6.1.1. ANALIZA SUROGAT PODATAKA\r\nImajući u vidu činjenicu da je priroda mehanizma nastanka rudarskih potresa (deterministička ili stohastička), čije se ponašanje izučava ispitivanjem odgovarajućih vremenskih serija, najvažniji preduslov primene adekvatne analize i validne interpretacije rezultata, tendencija savremenih istraživanja je da analizi vremenskih serija gotovo uvek prethodi analiza surogat podataka, što je izvedeno i u okviru ovog istraživanja.\r\nGrubo posmatrano, analizom surogat podataka vrši se određena transformacija posmatrane vremenske serije, tako da ona približno odgovara nekoj od poznatih funkcija, pa se potom sprovodi niz statističkih testova, kojima se postavljene nulte hipoteze o mogućoj linearnoj ili nelinearnoj prirodi vremenskih serija prihvataju ili odbacuju. Prema tome, nulte hipoteze predstavljaju početnu pretpostavku da je posmatrana vremenska serija ili potpuno nasumična, ili da pripada određenoj klasi linearnih procesa. Provera ovih nultih hipoteza operativno se sprovodi tako što se formira surogat vremenska serija, i to transformacijom početne vremenske serije (nasumičnom permutacijom, ili izmenama faze ili amplitude spektra snage), i upoređivanjem dobijenih surogat serija sa početnom vremenskom serijom. Kao statistički test, za upoređivanje početne i surogat vremenskih serija koristi se jednostavni algoritam nelinearnog predviđanja [Kantz, Schreiber, 2004].\r\nPrema pojedinim autorima, analiza surogat podataka predstavlja i deo ,,preprocessing-a'', odnosno pripada postupcima kojima se posmatrana vremenska serija ,,priprema'' za nelinearnu analizu. Naziva se još i tehnikom ,,reuzorkovanja'' ili ponovnog uzorkovanja, s obzirom na to da se odgovarajućim transformacijama, od početne vremenske serije formiraju nizovi surogat podataka, koji se potom upoređuju sa podacima registrovanja oscilovanja tla [Kugiumtzis, Tsimpiris, 2010].\r\nOpšti princip ove analize može se opisati na sledeći način. Najpre se pretpostavlja da vremenska serija vodi poreklo od određene klase linearnih procesa, odnosno postavlja se nulta hipoteza da izučavani niz podataka pripada nekom od poznatih linearnih sistema. Potom se za svaku nultu hipotezu generiše surogat vremenska serija i vrši niz statističkih testova radi poređenja surogat i početne vremenske serije. Ukoliko postoji statistički značajna razlika između surogat niza i početne vremenske serije, nulta hipoteza se odbacuje. Pri postavljanju nulte hipoteze obično se polazi od jednostavnih pretpostavki ka sve kompleksnijim i sofisticiranijim modelima, ukoliko rezultati statističke analize ne zadovoljavaju nultu hipotezu.\r\nSurogat analiza podataka je neophodna pri analizi vremenskih serija, s obzirom na to da jednostavni haotični sistemi mogu pokazivati ponašanje slično nasumičnom (proizvoljnom), i obratno [Sprott, 2003]. U opštem slučaju, pristupa se testiranju tri nute hipoteze [Perc et al., 2008]:\r\n(1) podaci su nezavisni slučajni brojevi sa određenom, ali nepoznatom raspodelom;\r\n(2) podaci vode poreklo od stacionarnog stohastičkog procesa sa Gausovom raspodelom stohastičkog dela;\r\n(3) podaci vode poreklo od stacionarnog Gausovog procesa koji je izmenjen\r\nnepoznatom nelinearnom funkcijom.\r\nTestiranje ovih hipoteza vrši se na dovoljnom broju surogat vremenskih serija, kako bi nulte hipoteze sa sigurnošću mogle biti prihvaćene ili odbačene. U opštem slučaju, 20 surogat vremenskih serija se smatra dovoljnim brojem, pri čemu, ukoliko je više od jednog testa sa različitim ishodom, nulta hipoteza se odbacuje [Kantz, Schreiber, 2004].\r\nGenerisanje surogat vremenskih serija vršeno je u okviru Matlab programskog paketa MATS, razvijenog u radu [Kugiumtzis, Tsimpiris, 2010], a njihovo testiranje na osnovu algoritma zeorpred.m razvijenog na Mekalester koledžu u Mineanapolisu (SAD) \r\n■ Formiranje nultih hipoteza i nizova surogat podataka.\r\nZa svaku od navedenih nultih hipoteza generiše niz od 20 surogat vremenskih serija, koje se potom upoređuju sa početnom vremenskom serijom, a potom vrši statistička analiza njihove sličnosti ili razlike.\r\n1. I nulta hipoteza - vremenska serija predstavlja permutovanu (,,ispreturanu'') početnu vremensku seriju.\r\nU ovom slučaju, od originalne (početne) vremenske serije, nizovi surogat podataka dobijaju se nasumičnim mešanjem postojećih podataka (bez ponavljanja).\r\n2. II nulta hipoteza - vremenska serija vodi poreklo od stacionarnog stohastičkog procesa sa Gausovom raspodelom stohastičkog dela.\r\nSurogat vremenske serije u okviru druge nulte hipoteze dobijaju se primenjujući postupak fazne randomizacije [Kantz, Schreiber, 2004]. Naime, ukoliko posmatrana\r\nvremenska serija vodi poreklo od stacionarnog linearnog stohastičkog procesa sa Gausovom raspodelom stohastičkog dela (ARMA procesi), onda se surogat podaci formiraju Furijeovim transformacijama početne vremenske serije, pri čemu se faza transformacije za svaku surogat vremensku seriju razlikuje, što se postiže primenom tehnike fazne randomizacije. Preciznije, najpre se određuje diskretna Furijeova transformacija1 sk originalne vremenske serije sn:\r\nPotom se vrši ,,randomizacija'' diskretne Furijeove transformacije sk, i to množenjem sa nasumično izabranim fazama:\r\ngde su nasumično izabrane vrednosti iz intervala [0,2n]. Potom se pristupa određivanju inverzne Furijeove transformacije2 od sk':\r\nčime se dobija surogat vremenska serija sn . Ovaj postupak se potom ponavlja za svaku sledeću surogat vremensku seriju, ali svaki put sa različitom fazom pk iz intervala [0,2n].\r\n3. III nulta hipoteza - vremenska serija vodi poreklo od stacionarnog Gausovog procesa koji je izmenjen nepoznatom nelinearnom funkcijom sn:\r\ngde xn predstavlja ARMA3 proces, sa prvim sabirkom koji predstavlja autoregresivni deo (AR) i opisuje unutrašnju dinamiku procesa, i sa drugim sabirkom, koji se naziva pokretnim prosekom4 (MA), gde su nn nezavisni brojevi sa Gausovom raspodelom.\r\nS obzirom na to da je funkcija sn u ovom slučaju nepoznata, postupak formiranja surogat vremenskih serija (AAFT postupak - postupak fazne randomizacije sa prilagođenom amplitudom) sastoji se u sledećem. Najpre se podaci početne vremenske serije skaliraju tako da odgovaraju Gausovoj raspodeli (efekat nepoznate nelinearne funkcije). Potom se primeni postupak fazne randomizacije, opisan u prethodnom koraku. Nakon primene inverzne Furijeove transformacije, dobijeni surogat podaci se na kraju opet skaliraju tako da odgovaraju raspodeli početne vremenske serije. Ovaj postupak se ponavlja za svaku surogat vremensku seriju, ali svaki put sa različitom fazom u postupku fazne randomizacije.\r\nMeđutim, s obzirom na to da se inverznom Furijeovom transformacijom ne dobija potpuno rekonstruisana vremenska serija, pre svega u pogledu spektra snage (koji je sa većom količinom šuma), vrši se i prilagođavanje amplitude Furijeovih transformacija na sledeći način (IAAFT postupak - iterativni postupak fazne randomizacije sa prilagođenom amplitudom). Neka su | sk | željene amplitude Furijeove transformacije početne vremenske serije. Izvršimo postupak fazne randomizacije, kada dobijenu surogat Furijeovu transformaciju sk' dalje menjamo poređenjem amplituda Furijeove transformacije početne i surogat vremenske serije:\r\nOvaj korak se ponavlja sve dok se količina šuma u spektru snage na svede na prihvatljivu meru. Nakon izvedene inverzne Furijeove transformacije za sk'', dobijena surogat vremenska serija poredi se sa početnom vremenskom serijom.\r\n■ Testiranje hipoteza.\r\nTestiranje predloženih nultih hipoteza po pravilu se zasniva na proceni greške u predviđanju (predviđanje nultog reda5). Naime, formira se 20 surogat vremenskih serija od početne vremenske serije, a potom se za svaku od njih određuje greška u predviđanju za određeni broj koraka (10 jedinica). Nulta hipoteza se prihvata ukoliko je greška u predviđanju za početnu vremensku seriju veća od greške u predviđanju za surogat vremenske serije, ili su one u korelaciji, u 95% ispitanih slučajeva (19 od 20). Drugim rečima, ukoliko ovaj uslov nije ispunjen za više od jedne surogat vremenske serije, onda se nulta hipoteza odbacuje.\r\nU opštem slučaju, neka je s1,....sN skalarna vremenska serija, T vremenski interval razvijanja i m dimenzija razvijanja. Tehnikom razvijanja formiramo niz razvijenih vektora S(m-1)T+1,...,SN. Da bismo predvideli vrednost vremenske serije za An vremenskih jedinica unapred, biramo parametar e, koji je reda veličine rezolucije merenja. Naime, da bi naše predviđanje bilo adekvatno i dovoljno precizno, ono ne može biti zasnovano samo na jednom najbližem stanju, pošto su međurastojanja između dva susedna podatka sa određenom greškom (šumom) reda veličine e. Stoga, sve tačke u faznom prostoru u radijusu e oko naše početne tačke moraju biti uzete u obzir kao dobri „prediktori\". Umesto proizvoljnog izbora bilo koje od tačaka iz navedenog radijusa, uzima se aritmetička sredina individualnih predviđanja vrednosti tačaka iz radijusa e. Drugim rečima, formira se okolina UE(sN) radijusa e oko tačke sN. Za sve tačke snCUE(sN), to jest, za sve tačke na manjem rastojanju od e u odnosu na sN, vrši se predviđanje sN+^n. Konačno predviđanje vrednosti na kraju predstavlja aritmetičku sredinu svih predviđanja baziranih na poiedinačnim tačakama:\r\ngde l UE(sN) l označava broj elemenata u okolini UE(sN). U slučaju da nema elemenata u okolini sN bližih od e, onda povećavamo e sve dok ne nađemo pogodno rastojanje sa dovoljnim brojem tačaka za predviđanje.\r\nDrugim rečima, za sledeću (buduću) vrednost jedne određene tačke uzima se srednja vrednost budućih vrednosti svih tačaka sn u okolini U izabrane tačke sN.\r\nModifikovani pristup, u odnosu na prethodno predstavljen, podrazumeva beskonačno veliku okolinu U, ali sa težinskim rastojanjima susednih vrednosti:\r\ngde ^(l sN-sn l) predstavlja težinsko rastojanje tačaka sn od tačke za koju se vrši predviđanje sN.\r\nDa bismo ocenili grešku predviđanja, za svaki predviđeni vremenski korak An, pristupamo izračunavanju srednje kvadratne greške:\r\nkoja određuje razliku između predviđenih i stvarno merenih vrednosti.\r\n■ Verifikacija predloženog postupka testiranja hipoteza.\r\nU cilju verifikacije tačnosti predloženog postupka, koji podrazumeva određivanje greške u predviđanju početne (s) i surogat vremenskih serija (s0), izvršena je provera koristeći opšte poznate primere izrazito linearnog stohastičkog sistema (ljubičasti šum6) i nelinearnog sistema (Lorencov sistem7). Kao što je opšte poznato, ljubičasti šum predstavlja primer linearnog stacionarnog stohastičkog procesa, sa frekventnim spektrom takvim da je spektar snage obrnuto proporcionalna frekvenciji (stoga se često ovakav obojeni šum naziva i 1/f šum). Za razliku od belog šuma, ljubičasti šum predstavlja niz korelisanih vrednosti, koje nisu međusobno nezavisne. S druge strane, Lorencov sistem reprezentuje model konvektivnih strujanja vazduha u atmosferi, predstavljen sistemom običnih diferencijalnih jednačina, čija su rešenja ,,haotična'' za određene vrednosti parametara sistema.\r\nNizovi podataka koji simuliraju rešenja Lorencovog sistema diferencijalnih jednačina dobijeni su numeričkom integracijom u okviru programskog paketa WINPP. U cilju ispitivanja (ne)linearnosti, pristupilo se generisanju 20 surogat nizova podataka za svaki od ispitivanih sistema, a potom je izvršeno testiranje sve tri nulte hipoteze. Pre početka testiranja, za originalne (početne) nizove podataka (Lorencov sistem i ljubičasti šum), kao i za sve surogat nizove, bilo je potrebno odrediti optimalne vrednosti vremenskog intervala razvijanja, kao i minimalne vrednosti dimenzije razvijanja. Dok je postupak uzajamne informacije8 omogućio dobijanje optimalnog vremenskog intervala razvijanja, postupak ,,prividno najbliže (susedne) vrednosti''5 je dao rezultate jedino u slučaju ispitivanja originalnog Lorencovog sistema (m=3). U slučaju originalnog niza podataka koji predstavljaju ljubičasti šum, kao i u slučaju svih nizova surogat podataka za oba sistema, minimalna dimenzija razvijanja nije mogla biti određena, s obzirom na to da su rezultati analize ukazivali na povećanje procenta prividno najbližih vrednosti sa povećanjem dimenzije razvijanja. Međutim, probanjem se ustanovilo da nema kvalitativnih razlika u distribuciji greške predviđanja za vrednosti dimenzije razvijanja m=1-10, odnosno da je predloženi postupak određivanja greške u predviđanju robustan u odnosu na izbor minimalne dimenzije razvijanja, tako da je za sve analize usvojena vrednost m=1.\r\nKao prvi korak, izvedeno je testiranje nulte hipoteze (da je vremenska serija sastavljena od nasumično izabranih brojeva, bez međusobne korelacije). Rezultati testiranja prikazani na slici 6-1 omogućuju odbacivanje nulte hipoteze u oba posmatrana slučaja.\r\nSlika 6-1. Testiranje I nulte hipoteze: (a) ljubičasti šum; (b) Lorencov sistem. Očigledno je da u oba slučaja nulta hipoteza može biti odbačena, s obzirom na to da je e00,05), odnosno statistički beznačajna, a za koeficijent b - što manja (p^-0).\r\nTabela 6-1. Test korelacije distribucije greške predviđanja za 20 nizova surogat podataka i početnu vremensku seriju pri testiranju II hipoteze, kod ljubičastog šuma.\r\nKao što se iz Tabele 6-1 može videti, u 2 od posmatranih 20 slučajeva odsečak na y- osi (a) nije statistički zanemarljiv (p<0,05), tako da možemo da odbacimo nultu hipotezu. Rezultati testiranja treće hipoteze (podaci vode poreklo od stacionarnog Gausovog procesa koji je izmenjen nepoznatom nelinearnom funkcijom) prikazani su na slici 6-3. Kao i kod testiranja prethodne dve hipoteze za Lorencov sistem, nulta hipoteza se može odbaciti za surogat nizove podataka dobijene primenom ,,običnog'' i iterativnog postupka ispitivanja surogat nizova podataka sa ,,podešenom'' amplitudom (AAFT i IAAFT). S druge strane, nulta hipoteza ne može biti odbačena u slučaju ljubičastog šuma.\r\nSlika 6-3. Testiranje III hipoteze; (a) ,,običnim'' i (b) iterativnim postupkom ispitivanja nizova surogat podataka sa ,,podešenom'' amplitudom (AAFT i IAAFT). Leva kolona predstavlja rezultate ispitivanja ljubičastog šuma, desna kolona - Lorencovog sistema. Crvena linija označava grešku u predviđanju za početnu vremensku seriju (e0), dok crne linije označavaju distribuciju greške za nizove surogat podataka (e).\r\nKao i u slučaju testiranja II hipoteze, rezultati testiranja III hipoteze mogu se uporediti sa početnom vremenskom serijom pomoću testa korelacije, radi utvrđivanja statističke značajnosti međusobnog odstupanja. U Tabeli 6-2 date su vrednosti korelacionog koeficijenta (R2) i praga značajnosti (p) za uporedni test svakog od 20 surogat nizova (dobijenih AAFT i IAAFT postupkom) sa početnim nizom podataka. Rezultati testiranja III hipoteze AAFT postupkom pokazuju da se samo u jednom od posmatranih 20 slučajeva hipoteza H0 može odbaciti, s obzirom na to da odsečak na y- osi (a) nije statistički zanemarljiv (p<0,05). U ostalim slučajevima, prag značajnosti p>0,05, sa visokom vrednošću koeficijenta korelacije (R2=0,994-1), tako da se nulta hipoteza H0 ne može odbaciti. Prema tome, može se zaključiti da ljubičasti šum pripada klasi stacionarnih Gausovih stohastičkih procesa izmenjenih nepoznatom nelinearnom funkcijom. S druge strane, rezultati testiranja III hipoteze IAAFT postupkom pokazuju da se u tri od posmatranih 20 slučajeva hipoteza H0 može odbaciti, s obzirom na to da je prag značajnosti p<0,05. U ostalim slučajevima, prag značajnosti p>0,05, sa visokom vrednošću koeficijenta korelacije (R2=0,994-0,998), tako da se nulta hipoteza H0 može odbaciti.\r\n6.1.1. ANALIZA NELINEARNIH VREMENSKIH SERIJA\r\nMatematičko modelovanje prirodnih procesa i njihovih pojava najčešće zahteva egzaktan matematički model izučavanog sistema, aproksimiran do određenog stepena. Međutim, u inženjerskoj praksi, jako je čest slučaj da ne raspolažemo matematičkim modelima posmatranih procesa, čime se u velikoj meri otežava izučavanje dinamike sistema. Razlog tome leži u činjenici da sam mehanizam procesa ili nije dovoljno dobro izučen da bi bio precizno definisan, ili je mehanizam nastanka procesa toliko kompleksan da ga nije moguće modelovati relativno jednostavnim matematičkim modelima (običnim diferencijalnim jednačinama). U tom slučaju, potrebno je okrenuti se raspoloživim podacima, koji su, u inženjerstvu, najčešće predstavljeni u obliku velikog broja podataka osmatranja samog procesa ili njegovih efekata. Drugim rečima, potrebno je na osnovu podataka osmatranja suditi o dinamičkim svojstvima posmatranog procesa.\r\nTabela 6-2. Test korelacije distribucije greške predviđanja za 20 nizova surogat podataka i početnu vremensku seriju pri testiranju III hipoteze (AAFT i IAAFT), kod ljubičastog šuma.\r\nPrimenjene metode uglavnom su bazirane na analizi nelinearnih vremenskih serija ubrzanja oscilovanja tla za vreme potresa, uz potvrdu analizom surogat podataka. Analiza je izvršena koristeći teoremu razvijanja u određenom vremenskom intervalu10.\r\nIdeja o rekonstrukciji faznog prostora iz vremenske serije potiče od H. Vitnija11 (1936), koji je tvrdio da preslikavanje iz n-dimenzionalne mnogostrukosti u 2n+1- dimenzionalni Euklidski prostor predstavlja potapanje: ,,slika'' n-dimenzionalne mnogostrukosti je u potpunosti razvijena u većem prostoru. Preciznije rečeno, bilo koje dve tačke u n-dimenzionalnoj mnogostrukosti se ne mogu preslikati u istu tačku u 2n+1-dimenzionalnom prostoru. Ako 2n+1 nezavisnih signala sistema posmatramo kao preslikavanje skupa stanja u 2n+1 dimenzionalnom prostoru, Vitnijeva teroma o razvijanju implicira da svako stanje može jednoznačno biti određeno vektorom dužine 2n+1, čime se, u suštini, rekonstruiše fazni prostor. Potom su, gotovo istovremeno, ali međusobno nezavisno, grupa sa Univerziteta u Kaliforniji, Santa Kruz, predvođena Pakardom [1980] i Dejvid Ruel12 predložili ideju korišćenja ,,vremenski pomerenih'' koordinata za rekonstrukciju faznog prostora posmatranog dinamičkog sistema. Naime, oni su pošli od činjenice da, ukoliko želimo da nađemo prvi izvod promenljive s(n), oblika:\r\nonda moramo da izvršimo aproksimaciju konačnom razlikom, pošto smo merenja izvršili samo za vremenski interval TS. Aproksimacija prvog izvoda data je u sledećem obliku:\r\nImajući u vidu konačnu vrednost TS, prethodni izraz predstavlja samo grubi visokopropusni filter merenih podataka. Kada bismo hteli da procenimo vrednost drugog izvoda ds2(t)/dt2, dobili bismo samo grubu aproksimaciju, slično i sa trećim izvodom, itd. Ukoliko malo podrobnije pogledamo prethodni izraz za aproksimaciju izvoda, videćemo da se pri svakom koraku već poznatoj informaciji o izmerenoj vrednosti s(n) dodaju merenja u sledećim vremenskim intervalima, dobijenim množenjem vrednosti TS. Glavna ideja proistekla iz ovog razmatranja je da izvodi uopšte nisu potrebni da bi se formirao koordinatni sistem u kojem bi se ,,otkrio'' raspored trajektorija u faznom prostoru, već da vremenski ,,pomerene'' promenljive s(n+r)=s(to+(n+T)rs) mogu biti direktno korišćene, gde je T celobrojna vrednost koja treba da se odredi.\r\nZatim je F. Takens [1981] pokazao da se rekonstrukcija može izvesti samo sa jednom merenom veličinom. Naime, Takens je dokazao da je, umesto generisanja 2n+1 signala, razvijanje signala u određenom vremenskom intervalu T, oblika: y(t)y(t-r)y(t-2T),...y(t-2m), dovoljno za razvoj n-dimenzionalne mnogostrukosti. Glavnu posledicu Takensove teoreme predstavlja mogućnost razvoja skalarne vremenske serije u višedimenzionalni fazni prostor. Pored ove teoreme, postoje i druge metode razvijanja skalarnih vremenskih serija u fazni prostor [Landa, Rosenblum, 1989; Mindlin i dr., 1991]. Međutim, razvijanje u određenom vremenskom intervalu je verovatno jedina sistematska metoda koja polazi od skalarnih podataka i vodi do multidimenzionalnog faznog prostora, i sigurno je najviše citirana i primenjivana u recentnim publikacijama. Sva testiranja u okviru analize nelinearnih vremenskih serija izvršena su u programskom paketu prof. M. Perca sa Univerziteta u Mariboru [Perc, 2012].\r\nIzvedena analiza sastojala se iz tri sukcesivna koraka:\r\n• I korak - analiza surogat podataka;\r\n• II korak - izvođenje preliminarnih testova (konstrukcija faznog portreta i Furijeovog spektra snage);\r\n• III korak - primena teoreme o razvijanju skalarne vremenske serije u faznom prostoru.\r\n■ PRELIMINARNI TESTOVI\r\nPre početka sprovođenja analize, vrlo često je korisno izvesti nekoliko jednostavnih testova, koji mogu kvalitativno ukazivati na dinamiku posmatranog sistema. U opštem slučaju, primenjuju se dva osnovna postupka:\r\na. konstrukcija faznog portreta sistema - u slučaju rudarski generisanih potresa, na apscisu nanosimo vrednosti pomeranja, a na ordinatu vrednosti brzine\r\noscilovanja tla za vreme potresa. Na osnovu rasporeda trajektorija (rešenja) sistema, moguće je dati kvalitativnu ocenu dinamike sistema; b. izračunavanje Furijeovog spektra snage - Furijeov spektar snage predstavlja vrlo koristan postupak kod ocene jednostavnih dinamičkih sistema:\r\n• ravnotežno stanje predstavljeno je neznatnom fluktuacijom amplitude, koja već za male vrednosti frekvencija konvergira ka nultoj vrednosti;\r\n• periodično (oscilatorno) ponašanje pokazuje jedan ,,pik'' u spektru snage (s obzirom na to da ima jednu izraženu frekvenciju);\r\n• kvaziperiodično ponašanje (torus) pokazuje najmanje dva jasno izražena ,,pika'' u spektru snage, bez ikakvog kontinualnog šuma između pikova. U slučaju postojanja kvaziperiodičnog ponašanja u više dimenzija, broj pikova može biti veći od 2;\r\n• deterministički haotično ponašanje se manifestuje širokopojasnim, kontinulanim šumom u spektru snage, bez jasno izraženih ,,pikova''.\r\n■ RAZVIJANJE VREMENSKE SERIJE U FAZNOM PROSTORU\r\nNakon primene preliminarnog testova pristupa se analizi nelinearnih vremenskih serija. Ključ ove analize leži u razvijanju skalarne vremenske serije u višedimenzionalnom faznom prostoru. Međutim, s obzirom na to da se, u ovom slučaju, osmatra samo jedna promenljiva (ubrzanje tla u vertikalnom pravcu z za vreme potresa), postavlja se pitanje validnosti analize. Drugim rečima, sistem nije u potpunosti definisan, tako da ni rezultati analize ne mogu biti pouzdani. Tako, u slučaju rudarski generisanih potresa, oscilovanje tla je u potpunosti definisano sa tri promenljive - ubrzanje u verikalnom pravcu i dva međusobno upravna horizontalna pravca (sistem sa tri stepena slobode). Imajući u vidu prethodno, analiza polazi od pretpostavke da su sve promenljive determinističkog dinamičkog sistema suštinski povezane, odnosno utiču jedna na drugu (međusobno spregnut sistem), što zapravo predstavlja suštinu Takensove teoreme o razvijanju. Posledica ove tvrdnje je sledeća: ako je u vremenskom trenutku t samo vrednost promenljive z poznata, onda sledeća izmerena vrednost z u vremenskom trenutku t+T mora implicitno sadržati određenu količinu informacije i o promenljivama x i y. Nastavljajući merenje promenljive z u vremenskim trenucima t+2R, t+3R,.., istovremeno prikupljamo informacije i o drugim dvema promenljivama, x i y. Zapravo, ukoliko vremenski interval T pravilno odaberemo, onda su podaci o vrednostima promenljivih x i y dovoljni, i možemo vrednosti promenljive x u vremenskim trenucima t+T i t+2T koristiti kao njihovu zamenu. Naš konačni cilj je da, primenom ove tehnike, izvršimo rekonstruciju faznog prostora, i to razvijanjem postojeće skalarne vremenske serije z u vremenskom intervalu t+2T. Da bismo to uradili, neophodno je da najpre odredimo minimalnu dimenziju i optimalni vremenski interval razvijanja13. Minimalna dimenzija razvijanja određuje se najčešće tehnikom ,,prividno najbliže (susedne) vrednosti''14, kada dve početno bliske vrednosti promenljive moraju ostati bliske nakon izvršenog razvijanja. Drugim rečima, dve udaljene vrednosti promenljive pre razvijanja ne smeju postati bliske nakon izvršene rekonstrukcije, jer bi to ukazivalo na činjenicu da skalarna vremenska serija nije u potpunosti razvijena. Shodno tome, ova metoda omogućuje određivanje minimalne dimenzije razvijanja, koja predstavlja onu vrednost za koju je procenat prividno najbližih vrednosti blizak nuli. S druge strane, optimalni vremenski interval razvijanja određuje se primenom metode ,,uzajamne informacije''15. Naime, uzajamna informacija kvantifikuje količinu podataka o promenljivoj xt+T pod uslovom da je poznata vrednost xt. Smatra se da prvi lokalni minimum funkcije uzajamne informacije predstavlja optimalni vremenski interval razvijanja, imajući u vidu činjenicu da u tom slučaju xt+T daje najveću količinu informacija onome što već znamo na osnovu vrednosti promenljive xt, bez potpunog gubitka korelacije između te dve vrednosti. Izračunavanjem minimalne dimenzije i optimalnog vremenskog intervala razvijanja, sada smo u mogućnosti da rekonstruišemo razvijeni fazni prostor skalarne serije z. U tom prostoru dalje možemo sprovesti deterministički test (radi ispitivanja determinističnosti izučavanog procesa) i test stacionarnosti (radi ispitivanja stacionarnosti procesa), koji su potrebni uslovi da bi izučavani sistem pokazivao deterministički haotično ponašanje. Takođe, možemo pristupiti izračunavanju maksimalnog Ljapunovljevog eksponenta, kao mere eksponencijalne divergencije dve početno bliske susedne vrednosti.\r\nPrema tome, osnovne faze pri razvijanju skalarne vremenske serije u faznom prostoru su:\r\n• određivanje optimalnog vremenskog intervala razvijanja;\r\n• određivanje minimalne dimenzije razvijanja;\r\n• rekonstrukcija razvijenog faznog prostora;\r\n• test determinističnosti;\r\n• test stacionarnosti;\r\n• određivanje vrednosti najvećeg Ljapunovljevog eksponenta.\r\n■ ODREĐIVANJE OPTIMALNOG VREMENSKOG INTERVALA RAZVIJANJA\r\nPrvi korak analize vremenskih serija rudarski generisanih potresa predstavlja određivanje optimalne vrednosti vremenskog intervala razvijanja. Kao što je već prethodno rečeno, neophodno je odrediti vremenski interval između dve susedne vrednosti posmatrane promenljive, koji je dovoljno veliki, a da pri tom ne dođe do potpunog gubitka međusobne korelacije posmatranih vrednosti. Ukoliko je T isuviše malo, onda će dve susedne vrednosti biti toliko blizu da ih ne možemo razlikovati. Slično tome, ukolilo je T suviše veliko, onda su dve susedne vrednosti, sa statističkog stanovišta, potpuno nezavisne jedna od druge [Perc, 2006]. Prema [Fraser, Swinney, 1986], za beskonačan niz podataka, bez šuma, vrednost T može biti proizvoljno izabrana - drugim rečima, izabrana vrednost vremenskog intervala bitno ne utiče na dalji tok analize. Međutim, u praksi, vremenske serije su uvek ograničene, i sa određenom količinom šuma, koji ,,maskira'' koristan signal. Stoga je neophodno usvojiti određeni postupak koji omogućava utvrđivanje optimalne vrednosti parametra T. U opštem slučaju, primenjuju se dva postupka - konstrukcija autokorelacione funkcije i tehnika ,,obostrane (uzajamne) informacije''.\r\nAutokorelaciona funkcija, u suštini, meri korelaciju signala sa samim sobom, pomerenim za neki vremenski interval T, i definisana je na sledeći na čin:\r\ngde je t vremenska jedinica (korak), a T ukupni vremenski interval. Smatra se da ona vrednost vremenskog intervala za koju autokorelaciona funkcija prvi put uzima negativnu vrednost, ili opada kao 1/e, predstavlja optimalnu vrednost. Ovaj postupak se pokazao dobrim za regularne vremenske serije, dok za izrazito nelinearne serije, primena autokorelacione funkcije može dovesti do nepouzdanih rezultata, s obzirom na to da se bazira isključivo na linearnoj statistici, i ne uzima u obzir nelinearne korelacije.\r\nDrugi postupak, koji se i primenjuje u ovom radu, jeste tehnika ,,obostrane (uzajamne) informacije'' originalno predložena u radu [Fraser, Swinney, 1986], u kojem su autori predložili prvi minimum funkcije uzajamne informacije između vrednosti promenljive xt i xt+T kao optimalnu vrednost vremenskog intervala razvijanja. Naime, uzajamna informacija između vrednosti promenljive xt i xt+T određuje količinu informacije o stanju promenljive u trenutku t+T, pod uslovom da je stanje promenljive poznato u trenutku t. U opštem slučaju, pretpostavimo da imamo dva sistema, A i B, čije ponašanje osmatramo u nekom vremenskom intervalu. Neka su at i bk rezultati merenja prvog odnosno drugog sistema. Količina informacija (u bitovima) o merenju at na osnovu merenja bk data je, na osnovu informacione teorije [Gallager, 1968]:\r\ngde je PA(a) verovatnoća osmatranja vrednosti a u skupu svih vrednosti sistema A, PB(b) - verovatnoća osmatranja vrednosti b u skupu svih vrednosti sistema B, i P^fab) zajednička verovatnoća merenja a i b. Veličina I naziva se uzajamna informacija dva merenja at i bk, i simetrična je, u smislu količine informacija koju možemo saznati o bk na osnovu merenja at. Prosečna uzajamna informacija između merenja bilo koje vrednosti at sistema A i vrednosti bk sistema B predstavlja srednju vrednost svih mogućih merenja, lAB(aubk):\r\nU kontekstu posmatrane vremenske serije, s(n), neka početna vremenska serija s(n) predstavlja skup A, a rezultati merenja u vremenskom trenutku t+T, s(n+R) - skup B. Prosečna uzajmna informacija između osmatranja u trenutku n i n+T, odnosno prosečna količina informacija o s(n+R) na osnovu osmatranja s(n) može se izraziti na sledeći način:\r\nKada su merenja sistema A i B potpuno nezavisna, tada je \r\nNakon određivanja distribucije I(T) postavlja se pitanje: koju vrednost I(T) možemo koristiti za rekonstrukciju skalarne vremenske serije u faznom prostoru? U radu [Fraser, Swinney, 1986] predlaže se izbor one vrednosti TM za koju funkcija I(T) dostiže svoj prvi minimum. Vremenski interval TM izabran je kao interval u kojem su merene vrednosti skoro nezavisne, ali ne i statistički nezavisne:\r\nPrethodno navedeni kriterijum je intuitivnog i empirijskog karaktera, predložen na osnovu numeričkih kalkulacija većeg broja nelinearnih dinamičkih sistema. Kasnije, u radu [Liebert, Schuster, 1989] potvrđen je kriterijum ,,prvog minimuma'', putem izračunavanja korelacionog integrala, koji pokazuje zadovoljavajuće numeričke rezultate za vrednost TM.\r\n■ MINIMALNO PRIHVATLJIVA DIMENZIJA RAZVIJANJA\r\nSledeći korak analize vremenskih serija predstavlja određivanje optimalne vrednosti dimenzije razvijanja, s obzirom na to da je primarni cilj teoreme razvijanja da ,,obezbedi'' dovoljno veliki Euklidski prostor Rd, tako da niz tačaka dimenzije dA može biti potpuno razvijen. Drugim rečima, dve početno bliske vrednosti posmatranog niza tačaka (vremenske serije) moraju ostati bliske (susedne) i nakon razvijanja u faznom prostoru. Naime, niz vrednosti osmatranja ispitivanog sistema možemo smatrati trajektorijama (rešenjima) sistema. S obzirom na to da trajektorije na atraktoru ne smeju da se presecaju, to znači da posmatranu vremensku seriju moramo predstaviti u dovoljno velikom faznom prostoru u kojem će taj uslov biti ostvaren. Na primer, ukoliko je naš atraktor trodimenzionalan, može se desiti da u dve dimenzije i dalje postoje delovi prostora u kojem trajektorija sistema preseca samu sebe. Tako da ako povećamo dimenziju na 3, ispunićemo uslov o jedinstvenosti rešenja.\r\nNeophodno je naglasiti da će uslov o jedinstvenosti rešenja (,,nepresecanju'') biti ispunjen i za sve veće dimenzije faznog prostora od prve koja ispunjava ovaj uslov. Sa čisto teorijskog stanovišta, ovo je potpuno prihvatljivo, tako da će bilo koja dimenzija veća od minimalno prihvatljive poslužiti za dalja izračunavanja. Međutim, sa aspekta primene dobijenih vrednosti, rad sa dimenzijom faznog prostora većom od potrebne može uzrokovati dve vrste problema: (1) vreme pojedinih izračunavanja povećava se eksponencijalno sa povećanjem dimenzije razvijanja; (2) u prisustvu ,,šuma'' (nekorisnog signala), sav nepopunjen deo faznog prostora (pošto je njegov manji deo već u potpunosti pokriven atraktorom) biće ispunjen šumom i drugim tipovima nekorisnog signala. Imajući u vidu prethodno, pored određivanja potrebne dimenzije razvijanja, potrebno je odrediti i njenu minimalno dozvoljenu vrednost.\r\nTakensova teroema predlaže da je prihvatljiva dimenzija razvijanja: dE>2dA, gde je dE optimalna dimenzija razvijanja, a dA dimenzija atraktora. Ovaj kriterijum uvek omogućava dobijanje optimalne vrednosti, ali ne i najmanje prihvatljive, a što je jedan od dva uslova koji mora biti ispunjen.\r\nU opštem slučaju, postoje četiri metode za određivanje potrebne dimenzije razvijanja: (1) dekompozicija kovarijantne matrice16; (2) ,,zasićenje'' (saturacija) invarijanti sistema17; (3) metod prividno najbližih vrednosti18, i (4) metod ,,pravih'' vektorskih polja19. U daljem radu, biće primenjen postupak prividno najbližih vrednosti [Kennel i dr., 1992].\r\nOpšte posmatrano, postupak određivanja optimalne dimenzije razvijanja pomoću tehnike prividno najbliže vrednosti je prilično intuitivan - dve bliske vrednosti u posmatranoj vremenskoj seriji moraju ostati bliske nakon razvijanja; inače, u suprotnom, predložena dimenzija faznog prostora nije dovoljna za potpuno razvijanje atraktora. Naime, u dimenziji d, svaki vektor y(k)=[s(k), s(k+T), ...,s(k+(d-1)T)] ima najbliži vektor ^(k), pri čemu se pod blizinom podrazumeva neka funkcija rastojanja, poput Euklidskog rastojanja, koje se može izraziti na sledeći način:\r\nRd(k) je verovatno malo u slučaju velikog broja podataka, a za skup od N podataka, Euklidsko rastojanje je reda veličine 1/N^1^. Za dimenziju d+1, rastojanje između najbližih vektora je:\r\nPotrebna je neka granična vrednost RT kao kriterijum za određivanje procenta prividno najbližih vrednosti. Imajući to u vidu, ukoliko je:\r\nonda su najbliže vrednosti u trenutku k (ili to+kTs) prividne ili lažne. U praksi, za vrednosti 10< R T < 50, broj prividno najbližih vrednosti je konstantan.\r\nMeđutim, ovaj kriterijum ima i jedan vrlo važan nedostatak. Problem se javlja zbog činjenice da, radeći sa konačnim skupom podataka, sve tačke če na kraju imati svoje ,,najbliže tačke'', čije se rastojanje ne menja značajnije sa povećanjem dimenzije razvijanja. Međutim, kao što je već prethodno rečeno, to što su dve tačke u jednoj dimenziji bliske ne znači da su bliske i na samom atraktoru. Drugim rečima, bliskost dve tačke određuje dimenzija atraktora RA. Tako da, ukoliko dve tačke nisu na bliskom rastojanju, Rd(k)~RA, onda će rastojanje Rd+i(k) biti najmanje oko 2Rd(k), prema već pomenutuoj Takensovoj teoremi, što dovodi do drugog kriterijuma za određivanje procenta prividno najbližih vrednosti:\r\nU ovoj disertaciji, primenjuje se prvi kriterijum, zadavanjem vrednosti praga u predloženom intervalu 10< RT < 50.\r\n■ DETERMINISTIČKI TEST\r\nPo određivanju optimalnog vremenskog intervala i minimalne dimenzije razvijanja, pristupa se ispitivanju determinističnosti/stohastičnosti početnog sistema (ubrzanja oscilovanja tla u opštem slučaju) čije manifestacije izučavamo putem ispitivanja odgovarajućih vremenskih serija. Deterministički karakter početnog sistema, uz stacionarnost, predstavlja potreban preduslov da bi sistem pokazivao determinističkp haotično ponašanje. U ovoj disertaciji, za ispitivanje determinističnosti, korišćen je postupak determinističkog testiranja, originalno predložen u radu [Kaplan, Glass, 1992], a kasnije detaljnije obrazložen i proveravan na konkretnim primerima od strane M. Perca [2006].\r\nDeterministički test polazi od pretpostavke da je početni sistem, čije ponašanje izučavamo na osnovu vremenske serije promene određenog parametra, u osnovi deterministički, što implicira da se ponašanje sistema dalje može opisati sistemom običnih diferencijalnih jednačina prvog reda. Važna posledica ovog svojstva, koja sledi iz teorije o običnim diferencijalnim jednačinama, jeste da rešenja takvih sistema postoje i da su jedinstvena20. Jedinstvenost rešenja predstavlja važno svojstvo na kojem počiva deterministički test. Još jedna posledica ovog svojstva sistema je da, u slučaju postojanja diferencijalnih jednačina, njihova rešenja mogu biti predstavljena u obliku vektorskog polja. Međutim, u ovom slučaju, ne postoji definisan sistem jednačina koji opisuje mehanizam nastanka potresa (bar ne egzaktno)21. U tom slučaju, jedinstvenost rešenja sistema, čije se vremenske serije izučavaju, ne može biti direktno testirana. Ono što nam deterministički test omogućava jeste konstrukcija vektorskog polja rešenja sistema direktno iz vremenske serije, a potom i provera ispunjenosti uslova o jedinstvenosti rešenja.\r\nKonstrukcija vektorskog polja sastoji se u sledećem - rekonstruisan fazni prostor, zajedno sa trajektorijom (rešenjem) sistema, biva najpre podeljen na niz polja jednake veličine, istih dimenzija kao i minimalno prihvatljiva dimenzija razvijanja. Svakom polju kroz koji prolazi trajektorija sistema dodeljuje se jedan vektor na sledeći način: svaki prolaz i-te trajektorije kroz k-to polje generiše jedinični vektor et, čiji je pravac i smer određen pravcem i smerom trajektorije sistema unutar posmatranog polja. Konačna aproksimacija vektora Vk u k-tom polju faznog prostora predstavlja, u suštini, srednju vrednost svih vektora Vk generisanih svakim prolazom trajektorije kroz odabrano polje:\r\ngde je Pk broj svih prolaza trajektorije kroz k-to polje. Ponavljajući ovaj postupak za sva polja faznog prostora obuhvaćenih trajektorijom sistema, dobija se rekonstruisani pravac vektorskog polja. Ovde je reč ,,pravac'' naglašena, s obzirom na to da se ne raspolaže podacima o brzini prolaza trajektorije kroz određeno polje. Stoga, ne mogu se pretpostaviti ni apsolutne dužine pojedinačnih vektora, što nije ni potrebno za deterministički test. Ono što jeste važno je da, ukoliko vremenska serija vodi poreklo od determinističkog procesa, i ukoliko je podela faznog prostora na polja dovoljno precizna, svako rekonstruisano vektorsko polje sadržaće isključivo po jedan vektor. U ovom konkretnom slučaju, to znači da se vektori ne smeju presecati, pošto, u suprotnom, ne bi bio ispunjen uslov o jedinstvenosti rešenja. U slučaju presecanja vektora, to bi značilo da za jednu istu tačku sistema (koja je simbolički predstavljena poljem) postoje dva rešenja - u tom slučaju, za datu tačku unutar polja ne postoji jedinstveni vektor.\r\nRadi razvijanja statističkog kriterijuma za razlikovanje determinističke dinamike od dinamike koja je rezultat stohastičkog procesa, u radu [Kaplan, Glass, 1992] predloženo je poređenje sa slučajnim hodom22 u d dimenzija, predstavljenim veličinom Rkm, koja se\r\nizražava u zavisnosti od broja prolaza trajektorije kroz određeno k-to polje na sledeći način:\r\ngde je cm konstanta koja zavisi od minimalno prihvatljive dimenzije razvijanja i njena\r\n1 /2 1 /2 1 /2 1 /2 vrednost iznosi n , 4/(6^)1,z, 3n /321/z, za vrednost minimalne dimenzije razvijanja\r\nm=2, 3, i 4, redom.\r\nKao konačna mera determinizma, predložena je srednja vrednost Vk u odnosu slučajan hod (kao stohastički proces), nazvana determinističkim faktorom K, koji se izračunava na sledeći način:\r\ngde A predstavlja ukupan broj polja obuhvaćenih trajektorijom (rešenjem) sistema. Iz prethodnog izraza je jasno da ukoliko je VK=1, za potpuno deterministički sistem, onda je i vrednost K jedinična. Za VK<1, vrednost K više nije jedinična.\r\nMeđutim, treba biti svestan činjenice da je za realno osmatrane podatke K gotovo uvek manje od jedinice, na šta je ukazivano i u radu [Kaplan, Glass, 1992] smatrajući da je K=1 samo za beskonačni niz podataka. S obzirom na to da se u inženjerskoj praksi uvek raspolaže sa vremenski ograničenim skupom podataka, K nikad neće dostići jediničnu vrednost. U cilju donošenja ocene koja vrednost K bi mogla biti uzeta u obzir kao merodavna za ocenu determinističnosti sistema, sprovodi se niz testova, koji uglavnom podrazumevaju promenu vrednosti determinističkog faktora sa promenom optimalne dimenzije i vremenskog intervala razvijanja. Odgovarajući testovi biće primenjeni i u ovom slučaju, kod prikaza rezultata istraživanja, u okviru dela ,, Analiza vremenskih serija rudarski generisanih potresa''.\r\n■ TEST STACIONARNOSTI\r\nSledeći korak u analizi vremenskih serija predstavlja ispitivanje stacionarnosti početnog sistema, što, uz determinističnost, predstavlja osnovni preduslov\r\ndeterminističkog haotičnog ponašanja. Pod stacionarnošću sistema podrazumeva se svojstvo sistema da njegovi različiti delovi imaju istu zajedničku raspodelu verovatnoće bez obzira na vremenske ili prostorne transformacije. Drugim rečima, statistički pokazatelji poput srednje vrednosti i varijanse ne menjaju svoju vrednost pri promeni vremenskog intervala ili položaja. Slab stacionarni proces odlikuje se konstantnom srednjom vrednošću i varijansom, dok su kod pravog (jakog) stacionarnog procesa i svi momenti višeg reda jednaki (koeficijent asimetrije i mera spljoštenosti23) [Challis, Kitney, 1991].\r\nKod nelinearnih sistema, testiranje stacionarnosti po pravilu se zasniva na proceni greške u predviđanju. Grubo posmatrano, obično se usvajaju lokalna predviđanja, za razliku od globalnih, koji se najčešće koriste kod linearnih sistema. Postoje brojni predlozi različitih postupaka izučavanja lokalne dinamike [Pikovsky, 1986; Sugihara, May, 1990; Eckmann i dr., 1986; Farmer, Sidorowich, 1987]. Najjednostavniji i robustan pristup podrazumeva lokalnu aproksimaciju dinamike na osnovu određene aritmetičke sredine budućih vrednosti susednih ,,tačaka'', što je već opisano kod testiranja nultih hipoteza u okviru analize surogat podataka.\r\nKonkretno, prethodno opisani algoritam primenjuje se na sledeći način. Podelimo posmatranu vremensku seriju na niz segmenata jednake veličine, koji se ne preklapaju. Potom biramo dva segmenta S i Sj, pri čemu je S segment na osnovu koga se konstruiše model za predviđanje, a Sj neka bude test-segment - to jest, segment na kome se vrši testiranje formiranog modela. Potom se računa srednja kvadratna greška između predviđenih vrednosti segmenta Sj na osnovu segmenta Si3 i realno merenih vrednosti segmenta Sj. Ukoliko je greška veća od prosečne, onda vrednosti segmenta S očigledno predstavljaju ,,loš'' model za podatke u segmentu Sj. Drugo moguće objašnjenje tako velike greške može ležati u činjenici da je dinamika procesa koja se manifestovala segmentom S očigledno različita od dinamike procesa zabeležene segmentom Sj. I, kao treću mogućnost, može se izdvojiti nedovoljna gustina uzorkovanja u segmentu Si3 koji stoga ne daje dovoljnu količinu podataka za adekvatna predviđanja vrednosti u segmentu Sj. Naravno, najmanju (i gotovo nikakvu) grešku možemo očekivati za slučaj\r\ni=j, s obzirom na to da tada isti segment služi i za formiranje modela (treniranje) i za testiranje. Ovaj postupak se ponavlja za sve izdvojene segmente ispitivane vremenske serije.\r\n■ MAKSIMALNI LJAPUNOVLJEV EKSPONENT\r\nKao što je već opisano u poglavlju ,,Nelinearna dinamika i teorija haosa - pregled osnovnih pojmova'', maksimalni Ljapunovljev eksponent meri eksponencijalnu divergenciju između dve susedne početno bliske trajektorije sistema. Postoje dva postupka ocene vrednosti maksimalnog Ljapunovljevog eksponenta - Vulfova metoda i Rozenštajnova metoda. Obe metode se podjednako primenjuju u istraživanjima na temu nelinearne dinamike, tako da će njihove osnovne postavke biti objašnjenje u daljem tekstu.\r\nVulfova metoda\r\nVulfova metoda za određivanje veličine maksimalnog Ljapunovljevog eksponenta [Wolf i dr., 1985] u velikoj meri je intuitivna i direktna, i sastoji se u sledećem. Za odabranu tačku (vrednost) u faznom prostoru, p(t), nalazimo najbližu tačkup(i), koja je manja od nekog unapred usvojenog rastojanja s, \\p(i) - p(t)\\< s. Potom pratimo promenu rastojanja između dve odabrane tačke za neki fiksirani interval vremena v, koji bi trebalo da bude nekoliko puta veći od usvojene optimalne vrednosti vremenskog intervala razvijanja T. Ukoliko je sistem determinstički haotičan, onda je rastojanje dve tačke nakon vremena v, \\p(i+v) - p(t+v)\\ =sv, veće od početnog rastojanja s; u suprotnom, s ~ sv. Nakon vremenskog intervala v, traži se nova najbliža tačka, p(j), čije je rastojanje u odnosu na tačkup(t+v) jednako početnom s, pod uslovom da je ugaona razlika između vektora koji su obrazovani tačkama p(t+v) i p(i+v), odnosno p(t+v) i p(j) mala; u suprotnom bira se druga tačkap(j) koja ispunjava ovaj uslov. Postupak se ponavlja sve do poslednje tačke trajektorije sistema. Tada se maksimalni Ljapunovljev eksponent određuje na sledećin način:\r\ngde M predstavlja ukupan broj posmatranih vremenskih intervala. U ovom postupku, sukcesivne zamene nakon svakog usvojenog vremenskog intervala v ključne su za precizno određivanje maksimalnog Ljapunovljevog eksponenta. Ukoliko je fazni prostor pravilno konstruisan i potpuno razvijen, to znači da je dovoljno gusto popunjen trajektorijom sistema, tako da prethodno opisani algoritam daje pouzdane podatke. Međutim, ukoliko je trajektorija sistema slabo razvijena, bićemo primorani da prihvatimo dosta velike vrednosti početnih rastojanja, ili uglovnih razlika na kraju svakog vremenskog intervala v. U tom slučaju, ne može se tvrditi da je dobijena vrednost maksimalnog Ljapunovljevog eksponenta u potpunosti tačna.\r\nRozenštajnova metoda\r\nSlično Vulfovoj metodi, za procenu veličine maksimalnog Ljapunovljevog eksponenta, u radu [Rosenstein i dr., 1993] predlaže se sledeći postupak: najpre nalazimo sve bliske tačkep(k) koje su na manjoj udaljenosti od odabrane e u odnosu na izabranu referentnu tačku p(i). Na taj način dobijamo skup početnih tačaka za susedne trajektorije u faznom prostoru. U sledećem koraku, potrebno je izračunati prosečno rastojanje svih trajektorija u odnosu na referentnu trajektoriju, u funkciji određenog vremenskog intervala An, na sledeći način:\r\ngde s predstavlja broj različitih tačakap(k) koje ispunjavaju uslov da je \\p(k) -p(i) \\ 500) i različite vrednosti početnog rastojanja s. Drugim rečima, odabrano rastojanje s mora da bude što manje, ali i dalje dovoljno veliko tako da u proseku svaka referentna tačka ima bar nekoliko susednih (b>10).\r\nPrednost Rozenštajnove metode nad Vulfovom je u sledećem:\r\n- za razliku od Vulfove metode, Rozenštajnova metoda ne zahteva minimalno prihvatljivu dimenziju razvijanja kao ulazni parametar, što omogućava izučavanje i sistema za koje nije moguće odrediti dimenziju razvijanja;\r\n- procena veličine maksimalnog Ljapunovljevog eksponenta u ovom slučaju je intuitivna i direktna, pa samim tim su i dobijeni rezultati veće pouzdanosti.\r\nU ovoj disertaciji obe izložene metode primenjene su pri ispitivanju vremenskih serija ubrzanja oscilovanja tla.\r\n6.2. ANALIZA DINAMIČKOG PONAŠANJA FENOMENOLOŠKOG BARIDŽ- KNOPOF MODELA BLOKA SA OPRUGOM\r\nAnaliza dinamike BK modela izvedena je na tri načina: primenom bifurkacione analize, numerički (DDE-BIFTOOL) i posmatranjem ponašanja rešenja sistema. Bifurkaciona analiza obuhvata istraživanje kvalitativnih (topoloških) promena ponašanja dinamičkih sistema. Kao što je u poglavlju ,,Nelinearna dinamika i teorija haosa - pregled osnovnih pojmova'' već navedeno, bifurkacije predstavljaju one vrednosti kontrolnih parametara sistema, pri kojima dolazi do kvalitativne promene njihove dinamike. U opštem slučaju, promena dinamike sistema izučava se u uslovima variranja postojećih kontrolnih parametara sistema, ili uvođenjem novih parametara. Pri tome, izdvajaju se dva osnovna tipa bifurkacione analize:\r\n- lokalna bifurkaciona analiza, kojom se ispituju promene lokalne stabilnosti ravnotežnog položaja (fiksne tačke) variranjem kontrolnih parametara, kada sopstvene vrednosti rešenja ili realne vrednosti konjugovano-kompleksnih rešenja u fiksnoj tački menjaju predznak;\r\n- globalna bifurkaciona analiza, koja se primenjuje u slučaju kada se ,,veći'' dinamički skupovi, poput periodičnih orbita, ,,sudaraju'' sa fiksnom tačkom ili međusobno, pri čemu dolazi do kvalitativne promene dinamike sistema, koja se ne može osmatrati u blizini ravnotežnog položaja (fiksne tačke). Takve promene ne mogu biti razmatrane analizom stabilnosti u okolini fiksne tačke, odnosno ravnotežnog stanja.\r\nLokalna bifurkaciona analiza.\r\nU opštem slučaju, lokalna bifurkaciona analiza može biti izvedena analitičkim i numeričkim putem. U ovoj disertaciji, izvšena je analiza stabilnosti fenomenoloških modela nastanka rudarski generisanih potresa, putem izvođenja standardne lokalne bifurkacione analize (analitičkim putem), a potvrda dobijenih rešenja izvršena je numerički u okviru programskog paketa DDE-BIFTOOL. Bifurkaciona analiza je je primenjena pri izučavanju matematičkog modela R. Madariage, u radu [Erickson i dr., 2008], koji opisuje kretanje jednog bloka duž hrapave nepokretne podloge, a koji je oprugom određene krutosti povezan sa gornjom pločom koja uzrokuje kretanje sistema. Postojeći model je modifikovan uvođenjem vremenskog kašnjenja T u izraz za trenje.\r\nAnalitički postupak bifurkacione analize sastoji se u sledećem. U opštem slučaju, neka je x* fiksna (singularna) tačka sistema, i neka n(t)=x(t)-x* predstavlja mali otklon (poremećaj) od fiksne tačke x*. Da bismo analizirali ponašanje sistema u ovoj maloj okolini fiksne tačke, izvodimo diferencijalnu jednačinu za n po t:\r\njer je x* konstanta. Prema tome:\r\nRazvojem funkcije f(x) u Tejlorov red, dobijamo:\r\ngde O(n2) označava male kvadratne članove. Kako je f(x*)=0, jer je x* singularna tačka, sledi da je:\r\nKada jef'(x*)^0, kvadratni članovi O(n2) se mogu zanemariti, čime se dobija:\r\nKarakteristična jednačina sistema (6.29) može se predstaviti u opštem obliku:\r\ngde X predstavlja sopstvenu (karakterističnu) vrednost, od koje zavisi stabilnost sistema u maloj okolini fiksne (singularne) tačke. U opštem slučaju, X=oMfi, tako da se dva linearno nezavisna rešenja jednačine (6.29) mogu zapisati u sledećem obliku (primenom Ojlerove teoreme):\r\nU oba slučaja, realni deo sopstvene vrednosti X određuje stabilnost sistema, pri čemu možemo razlikovati nekoliko tipičnih slučajeva:\r\n- ukoliko obe sopstvene vrednosti X12 imaju realne delove sa negativnim predznakom, onda \\ q(t) \\ ^0 kada (fiksna tačka je stabilna);\r\n- ukoliko je bar jedna od sopstvenih vrednosti sa pozitivnim predznakom, onda postoji bar jedno rešenje n(t) takvo da \\ q(t) \\ ^^ kada (fiksna tačka je nestabilna-repeler, ili je sedlo);\r\n- ukoliko je za obe sopstvene vrednosti realni deo jednak nuli, odnosno imaju samo imaginarni deo, odgovarajuća rešenja osciluju u okolini fiksne tačke, odnosno niti joj se približavaju niti se udaljavaju od nje, za t^«. U tom slučaju, članovi višeg reda u jednačini (6.28) ne mogu se zanemariti.\r\nLinearizacija sistema se sprovodi zato što, prema Hartman-Grobmanovoj teoremi [Hartman, 1960; TpoGMaH, 1959] ponašanje dinamičkog sistema u okolini hiperboličke24 fiksne (ravnotežne) tačke je kvalitativno isto kao ponašanje linearizovanog sistema u blizini fiksne tačke, pod uslovom da sve sopstvene vrednosti linearizovanog sistema imaju realni deo Re(X)^0. Drugim rečima, izučavanje ponašanja dinamičkog sistema u maloj okolini ravnotežnog stanja (fiksne tačke) pogodnije je nakon izvršene linearizacije.\r\nNeophodno je naglastiti da se pre izvođenja linearizacije sistema oko fiksne tačke vrlo često pristupa bezdimenzionalizaciji posmatranog sistema jednačina, naročito ako je u pitanju sistem sa većim brojem parametara, čije bi razdvajanje u određene bezdimenzionalne grupe parametara u znatnoj meri olakšalo izvođenje postupka analize stabilnosti. Bezdimenzionalizacija matematičkog modela predstavlja konstruktivan način formulisanja modela samo pomoću bezdimenzionalnih veličina. Time se pre svega smanjuje broj parametara sistema. U opštem slučaju, ne postoji jedinstveni metod određivanja najpogodnjih bezdimenzionalnih veličina, već se za svaku jednačinu posebno konstruišu ove veličine putem pokušaja i pogreški25. U okviru poglavlja ,,Rezultati istraživanja'' biće prikazan postupak bezdimenzionalizacije, izveden za matematičke modele R. Madariage [2008] i Bekera [2000], koji opisuju kretanje jednog bloka.\r\nNakon izvedene linearizacije sistema, pristupa se određivanju karakteristične jednačine sistema i njenim rešavanjem po sopstvenoj (karakterističnoj) vrednosti X. Smenom dobijenih vrednosti X u karakterističnoj jednačini, određuju se relacije između posmatranih kontrolnih parametara26 sistema. Nakon toga, pristupa se crtanju bifurkacionih krivih, koje, u suštini, predstavljaju granice različitog dinamičkog ponašanja sistema u prostoru parametara. Crtanje bifurkacionih krivih je, u ovom slučaju, izvedeno korišćenjem programskog paketa MAPLE, i to u dvodimenzionalnom prostoru kontrolnih parametara. Dobijene bifurkacione krive proverene su korišćenjem usko specijalizovanog MATLAB paketa DDE-BIFTOOL (mr Igor Franović).\r\nProvera tačnosti određivanja bifurkacionih krivih vršena je simulacijom vremenskih serija i odgovarajućih faznih portreta za pojedine vrednosti kontrolnih parametara iznad i ispod određene bifurkacione krive.\r\nU cilju određivanja karaktera promene dinamike sistema pri prelazu preko bifurkacione krive, u narednom koraku se vrši izračunavanje izvoda po bifurkacionim krivama. Drugim rečima, u pojedinim odabranim tačkama duž krivih, dovoljne gustine,\r\nračunamo prvi izvod sopstvene vrednosti po odabranom kontrolnom parametru sistema. Dobijena pozitivna vrednost izvoda ukazuje na direktnu Hopfovu27 bifurkaciju, dok negativan izvod ukazuje na inverzni tip Hopfove bifurkacije. Naglašavamo da lokalnom analizom stabilnosti nismo u mogućnosti da odredimo tačan tip Hopfove bifurkacije (natkritična ili potkritična), s obzirom na to da globalno ponašanje sistema prevazilazi mogućnosti primenjene linearne analize.\r\nU petom koraku linearne analize stabilnosti, nakon linearizacije sistema, nalaženja karakteristične jednačine, pa nalaženja relacija između parametara sistema i izvoda po bifurkacionim krivama, pristupa se određivanju vrednosti kontrolih parametara sistema, za koje su rešenja sistema uvek stabilna. Ovo se može ostvariti primenom Rušeove teoreme, koju ove navodimo bez izvođenja dokaza.\r\nRušeova teorema: Ako su dve funkcije y(z) i p(z) analitičke unutar i na konturi C i \\ y(z) < y(z) | na C, onda y(z) i y(z)+ y(z) imaju isti broj rešenja unutar C.\r\nPri tome, najčešće se za funkciju $(z) usvaja netranscedentni deo karakteristične jednačine početnog sistema, dok se za funkciju y(z)+ p(z) uzima cela karakteristična jednačina sistema diferencijalnih jednačina sa kašnjenjem. Kao rezultat primene Rušeove teoreme dobija se odnos (oblast) kontrolnih parametara sistema, za koji su rešenja početne jednačine uvek stabilna.\r\nPoslednji korak u analizi stabilnosti predstavlja simulacija rešenja sistema za određene vrednosti kontrolnih parametara i za odgovarajući vremenski interval. Za potrebe izrade ove disertacije, simulacija rešenja je izvedena numerički, sastavljanjem posebnog programa u programskom jeziku C++ i njegovom implementacijom u Microsoft Visual Studio (mr Igor Franović).\r\nNumerička analiza ponašanja rešenja u parametarskom prostoru sistema\r\n■ Varijacije kontrolnih parametara.\r\nNumeričku analizu ponašnja rešenja u parametarskom prostoru sistema započinjemo najjednostavnijim postupkom - variranjem jednog ili više odabranih kontrolnih parametara, dok se vrednosti ostalih parametara drže konstantnim. U drugom koraku, pretpostavljamo da se promena vrednosti odabranih kontrolnih parametara može opisati nekom periodičnom funkcijom, opšteg oblika:\r\npri čemu amplitude tih oscilacija (perturbacija) moraju biti manje od vrednosti odabranog parametra r, Sr1m/s2, onda instrument nije registrovao sve komponente\r\nSlika 7-5. Broj potresa u funkciji magnitude, registrovanih u široj okolini rudnika ,,Rudna'' u periodu 1980-2002.g. Ukupan broj registrovanih potresa iznosi 15573 [Kwiatek, 2004]\r\n7.1.2. REZULTATI ANALIZE\r\nAnaliza ubrzanja oscilovanja tla izvršena je na primeru potresa registrovanog 2. februara 2001.g., na stanici ,,3 Maja'', instaliranoj u tlu. Energija potresa iznosila je 3,4x107J. Maksimalno registrovano ubrzanje tla (PGA) iznosilo je 51,6cm/s2, maksimalna brzina (PGV) - 0,649cm/s, a maksimalno pomeranje (PGD) 0,0442cm. Dužina trajanja potresa bila je tD=2,56s. Registrovanje oscilovanja tla vršeno je na svakih 0,0012s [Zembaty, 2004]. Za analizu su korišćene vremenske serije od 4050 podataka (pravac sever-jug), 3600 podataka (pravac istok-zapad) i 2750 podataka (vertikalni pravac). Na slici 7-6 prikazane su vremenske serije sve tri komponente pomeranja, brzine i ubrzanja oscilovanja tla.\r\nSlika 7-6. Vremenske serije oscilovanja tla (pomeranja, brzine i ubrzanja) za vreme potresa magnitude M~2 (3,4x107J) registrovanog 02.02.2001.g. u rudniku ,,Rudna'' u zapadnoj Poljskoj, na stanici ,,3 Maja'', u tri pravca: (x) - sever-jug; (y) - istok zapad i (z) - vertikalno [Zembaty, 2004].\r\nAnaliza je izvršena samo za komponente ubrzanja tla, s obzirom na to da maksimalno ubrzanje tla predstavlja jedan og glavnih parametara koji ulaze u proračun građevinskih konstrukcija [Boore, Bommer, 2005]. Takođe, seizmogram ubrzanja tla je mnogo ,,bogatiji'' u odnosu na registrovane promene brzine i, naročito, pomeranja tla. S druge strane, iako se ubrzanja u horizontalnom pravcu najčešće uzimaju u obzir, zbog njihove prirodne veze sa inercijalnim silama [Trifunac, Brady, 1975; Kramer, 1996], ovde će biti izvršena analiza sve tri komponente ubrzanja, radi ispitivanja dinamičkog ponašanja tla pri delovanju potresa u sva tri pravca. Pri tome, istovremeno sa analizom efekata rudarski izazvanih potresa, biće vršeno upoređivanje sa rezultatima analize tektonskih zemljotresa, a sve u cilju preciznijeg definisanja karaktera rudarski izazvanih potresa. Naime, za razliku od tektonskih zemljotresa, rudarski generisani potresi znatno kraće traju, u proseku do 5s. Shodno tome, manji je broj podataka za izvođenje analize, što može uticati na pouzdanost dobijenih rezultata. Dubine žarišta rudarski generisanih potresa su obično male, i kreću se u rasponu 2-4km. Takođe, vrlo brzo dolazi do prigušenja oscilacija u slučaju rudarski generisanih potresa. Tako na primer, ubrzanje oscilovanja tla od 0,39g izmereno na udaljenosti od 2,5km, redukuje se na 0,13g na rastojanju od 10km [Milford, Wium, 1991]. Na slici 7-7, na kojoj je dato poređenje registrovanog oscilovanja tla za vreme rudarski izazvanog potresa i tektonskog zemljotresa, pored dužine trajanja, jasna je i razlika u amplitudi signala.\r\nSlika 7-7. Poređenje registrovanog ubrzanja oscilovanja za vreme rudarski generisanog potresa (a), i Kraljevačkog zemljotresa (M5.4) od 03.10.2010.g. (b). Očigledna je razlika u dužini signala, amplitudi ubrzanja i karakteru oscilovanja.\r\n7.1.2.1. PRELIMINARNI TESTOVI\r\nPre početka izvođenja analize vremenskih serija, korisno je primeniti niz jednostavnih testova, u cilju dobijanja preliminarne slike dinamike izučavanog sistema. Prvi korak ka tome predstavlja konstrukcija faznog portreta za registrovani signal (slika 7-8). Kao što se može videti, fazni portret ne predviđa bogato dinamičko ponašanje. Naime, u slučaju rudarski generisanog potresa trajektorija dinamike sistema u faznom portretu ne pokazuje posebnu zakonitost ni u jednom pravcu oscilovanja tla, što može ukazivati na nasumično ponašanje2 ili deterministički haotično ponašanje.\r\nSlika 7-8. Fazni portret U-V (pomeranje-brzina) oscilovanja tla za vreme rudarski generisanog potresa, registrovanog na stanici ,,3 Maja'' u pravcu (a) sever-jug; (b) istok-zapad i (c) vertikalno. Očigledno je da se nikakva pravilnost u evoluciji trajektorija sistema ne može uočiti.\r\nSledeći korak preliminarne analize dinamičkog ponašanja podrazumeva izračunavanje Furijeovog spektra snage, kojim se naročito jasno može razdvojiti kvaziperiodično od haotičnog ponašanja [Sprott, 2003]. Slika 7-9 pokazuje primer Furijeovog spektra snage u jednom horizontalnom pravcu (sever-jug) i u vertikalnom pravcu. Kvalitativno sličan grafik se dobija i u drugom horizonalnom pravcu istok- zapad. Kao što se može videti sa slike 7-9, Furijeov spektar snage u oba predstavljena slučaja je kontinualan i dobro razvijen, bez izraženih pikova, što može ukazivati na postojanje haotičnog ponašanja. S druge strane, u radu [Sprott, 2003] ukazuje se na činjenicu da vrlo slični spektri snage mogu biti proizvod i potpuno nasumičnog ponašanja izučavanog sistema. S tim u vezi, sasvim je jasno zašto je potrebna gotovo ,,baterija'' testova kako bismo suzili moguće kandidate za ocenu dinamike oscilovanja tla tokom rudarski generisanog potresa.\r\nSlika 7-9. Furijeov spektar snage za ubrzanje tla registrovano na 3Maja stanici u (a) horizontalnom (sever-jug) i (b) vertikalnom pravcu. Kvalitativno sličan spektar snage se dobija i za ubrzanje tla u pravcu istok-zapad. Kontinualni spektar snage u oba slučaja ukazuje na moguće determinističkO haotično ili nasumično ponašanje posmatranog sistema.\r\n7.1.2.2. ANALIZA SUROGAT PODATAKA\r\nRezultati testiranja I hipoteze (da je vremenska serija sastavljena od nasumično izabranih brojeva, bez međusobne korelacije), gde su nizovi surogat podataka dobijeni\r\ngenerisanjem nasumično odabranih podataka početnog niza, prikazani su na slici 7-10. Jasno je da se u sva tri posmatrana slučaja nulta hipoteza H0 odbacuje, s obzirom na to da je e > e0 za svako n i za svaki ispitivani niz surogat podataka.\r\nSlika 7-10. Testiranje I nulte hipoteze: (a) pravac sever-jug; (b) pravac istok-zapad; (c) vertikalno. Crvena linija označava grešku u predviđanju za originalnu vremensku seriju (e0), dok crne linije označavaju distribuciju greške za nizove surogat podataka (e).\r\nNa slici 7-11 prikazani su rezultati testiranja druge hipoteze (podaci vode poreklo od stacionarnog stohastičkog procesa sa Gausovom raspodelom stohastičkog dela).\r\nSlika 7-11. Testiranje II nulte hipoteze: (a) pravac sever-jug; (b) pravac istok-zapad; (c) vertikalno. Crvena linija označava grešku u predviđanju za originalnu vremensku seriju (e0), dok crne linije označavaju distribuciju greške za nizove surogat podataka (e).\r\nI u ovom slučaju, pristupilo se generisanju 20 različitih nizova surogat podataka, primenom postupka fazne randomizacije. Za razliku od prethodnog slučaja, rezultati testiranja II hipoteze ukazuju na to da je e>e0 za pravac sever-jug i vertikalni pravac, dok je e0,05), odnosno statistički beznačajna, a za koeficijent b - što manja (p^-0).\r\nKao što se iz Tabele 7-2 može videti, u 3 od posmatranih 20 slučajeva za pravac sever-jug, odsečak na y-osi (a) je statistički zanemarljiv (p>0,05), sa visokom vrednošću koeficijenta korelacije (R2=0,99-1), tako da nultu hipotezu možemo odbaciti. Takođe, nultu ipotezu odbacujemo i za registrovana oscilovanja tla u vertikalnom pravcu, s obzirom na to da odsečak na y-osi ni u jednom slučaju nije beznačajan.\r\nPoslednji korak u analizi surogat podataka predstavlja testiranje hipoteze da vremenska serija vodi poreklo od stacionarnog Gausovog linearnog procesa, koji je modifikovan nepoznatom nelinearnom funkcijom. Kao u prethodnom slučaju, generisano je 20 surogata za svaki pravac oscilovanja tla, primenom ,,običnog'' i iterativnog postupka ispitivanja surogat nizova podataka sa ,,podešenom'' amplitudom (AAFT i IAAFT), i određena je greška predviđanja e (slika 7-12).\r\nTabela 7-2. Test linearne korelacije između distribucije greške predviđanja za 20 nizova surogat podataka i početnu vremensku seriju pri testiranju II hipoteze, za sve tri komponente ubrzanja.\r\nSlika 7-12. Testiranje III hipoteze; (a) pravac sever-jug; (b) pravac istok-zapad; (c) vertikalno. Leva kolona predstavlja rezultate ispitivanja ,,običnim''postupkom sa ,,podešenom'' amplitudom (AAFT), desna kolona - iterativnim postupkom ispitivanja nizova surogat podataka sa ,,podešenom'' amplitudom (AAFT i IAAFT). Crvena linija označava grešku u predviđanju za originalnu vremensku seriju (e0), dok crne linije označavaju distribuciju greške za surogat nizove podataka (e).\r\nRezultati testiranja III hipoteze ukazuju na to da jedino u slučaju oscilovanja tla (za slučaj testiranja IAAFT postupkom) u pravcu istok-zapad (slika 7-12b) ne možemo odbaciti nultu hipotezu (e<£o). U svim ostalim slučajevima, neophodno je izvođenje testa korelacije, na osnovu čega bi se pristupilo prihvatanju ili odbacivanju nulte hipoteze.\r\nU Tabeli 7-3 date su vrednosti korelacionog koeficijenta (R2) i praga značajnosti (p) za uporedni test svakog od 20 surogat nizova (dobijenih AAFT i IAAFT postupkom) sa početnim nizom podataka za pravac sever-jug i vertikalni pravac.\r\nU slučaju ispitivanja ubrzanja oscilovanja tla u pravcu sever-jug, rezultati testiranja III hipoteze AAFT postupkom pokazuju da je u pet od posmatranih 20 slučajeva odsečak na,y-osi (a) statistički zanemarljiv (p>0,05), sa visokom vrednošću koeficijenta korelacije (R2=0,99-1). Shodno tome, nulta hipoteza H0 može se odbaciti. S druge strane, rezultati testiranja III hipoteze IAAFT postupkom pokazuju da je u dva od posmatranih 20 slučajeva prag značajnosti p>0,05, sa koeficijentom korelacijem R2>0,998, tako da se nulta hipoteza H0 može odbaciti. Shodno tome, može se reći da se ubrzanje oscilovanja tla u pravcu sever-jug ne može svrstati u klasu vremenskih serija koje vode poreklo od stacionarnog Gausovog linearnog procesa, koji je modifikovan nepoznatom nelinearnom funkcijom.\r\nTabela 7-3. Test korelacije između distribucije greške predviđanja za 20 nizova surogat podataka i početnu vremensku seriju pri testiranju III hipoteze (AAFT i IAAFT), za registrovanja ubrzanja oscilovanja tla u pravcu sever-jug i u vertikalnom pravcu.\r\nU slučaju ispitivanja ubrzanja oscilovanja tla u vertikalnom pravcu, rezultati testiranja III hipoteze AAFT postupkom pokazuju da je u sedam od posmatranih 20 slučajeva odsečak na ,y-osi (a) statistički zanemarljiv (p>0,05), sa visokom vrednošću koeficijenta korelacije (R2=0,99-1), tako da se nulta hipoteza H0 može odbaciti. Takođe, rezultati testiranja III hipoteze IAAFT postupkom pokazuju da je u tri od posmatranih 20 slučajeva prag značajnosti p>0,05 za odsečak na y-osi (a), tako da se ubrzanje\r\noscilovanja tla u vertikalnom pravcu ne može svrstati ni u klasu vremenskih serija koje vode poreklo od stacionarnog Gausovog procesa, koji je modifikovan nepoznatom nelinearnom funkcijom.\r\nU slučaju ispitivanja ubrzanja oscilovanja tla u pravcu istok-zapad (Tabela 7-4), rezultati testiranja III hipoteze AAFT postupkom pokazuju da se u svim od posmatranih 20 slučajeva hipoteza H0 može odbaciti, s obzirom na to da odsečak na ,y-osi (a) nije statistički zanemarljiv (p<0,05), sa koeficijentom korelacije R2>0,973. Već je ranije naglašeno, da se nakon testiranja III hipoteze IAAFT postupkom za registrovano ubrzanje oscilovanja tla u pravcu istok.-zapad, hipoteza H0 ne može odbaciti, s obzirom na to da je e gubi za velike vrednosti argumenta i ^(0) = 1, dok uc predstavlja karakterističnu brzinu (Carlson and Langer, 1989). Radi pogodnosti, sistem (7.1.) može biti predstavljen u bezdimenzionalnom obliku, uvodeći nove promenljive:\r\nVeličina 2n/rnp predstavlja period oscilacija jednog bloka sa oprugom u odsustvu sile trenja klizanja (Carlson, Langer, 1989). Uvode nove promenljive (7.3), sistem (7.1) za jedan blok, u bezdimenzionalnom obliku, postaje (De Sousa Vieira, 1999):\r\ngde U predstavlja pomeranje bloka, v - bezdimenzionalnu brzinu bloka, a F predstavlja odgovarajući zakon trenja:\r\ngde je vc bezdimenzionalna karakteristična brzina. Kao što se može videti, sila trenja F između bloka i donje hrapave podloge zavisi samo od brzine bloka. U ovom slučaju, za razliku od pretpostavke u radu [Vieira, 1999], dozvoljeno je kretanje bloka ,,unazad''.\r\nPolazeći od jednačine (7.4) za model sa jednim blokom, uvodimo vremensko kašnjenje T i parametar jačine trenja c na sledeći način:\r\ngde su parametri isti kao u jednačini (7.1). Takođe, parametri T} i T2, odnosno cj and c2 uvode se i u model sa dva bloka (slika 7.24):\r\ngde k i k2 predstavljaju odnos sile u opruzi koja povezuje dva bloka, kc, i sila u oprugama koje povezuju blokove i gornju pokretnu ploču, kp, za prvi i drugi blok, redom. Parametri Vjc and predstavljaju bezdimenzionalne karakteristične brzine, dok su ostali parametri kao u jednačini (7.4).\r\nSlika 7-24. Baridž-Knopof model sa dva bloka (m predstavlja masu bloka, kc i kp su konstante krutosti opruga, V je brzina gornje pokretne ploče).\r\nDa bismo odredili fiksnu tačku sistema (ravnotežno stanje), ukoliko pođemo od\r\njednačine (7.7), i pretpostavimo da je Ui = U2 = v, Ui = U2 = 0 (De Sousa Vieira, 1999), a usvojimo da je k=k2=k, dobijamo:\r\nSabiranjem dve jednačine u sistemu (7.8) dobijamo:\r\na oduzimanjem druge od prve jednačine u sistemu (7.8) dobija se da je:\r\nodakle sledi:\r\nZamenjujući (7.11) u (7.9) dobijamo izraz za ravnotežno stanje sistema (7.7):\r\nU ovom slučaju, analiza je izvršena posmatranjem ponašanja rešenja u blizini ravnotežne tačke sistema. Pri tome, diferencijalne jednačine su rešavane numerički, primenom Runge-Kuta postupka četvrtog reda. Pri rešavanju diferencijalnih jednačina sa kašnjenjem, početna funkcija za U određena je u intervalu [-T, 0] jednačinom (7.6) za c=0 (za model sa jednim blokom) i jednačinom (7.4) za ci=c2=0 (za model sa dva bloka).\r\nRezultati su dobijeni variranjem kontrolnih parametara za model sa jednim blokom (c, T) i dva bloka (ci, c2, TI, i T2). Za svaki ispitani slučaj, prikazane su odgovarajuće vremenske serije i fazni portreti, zajedno sa određivanjem Furijeovog spektra snage i maksimalnog Ljapunoljevog eksponenta.\r\n7.2.1.1. MODEL SA JEDNIM BLOKOM\r\nU prvoj fazi analize, razmatra se početni model sa jednim blokom (7.4), za koji je usvojeno da je v=0,1 [Vieira, 1999]. Rezultati analize ukazuju na to da za različite vrednosti parametra vc, ovaj model pokazuje samo oscilatorno ponašanje različitih amplituda. Na slici 7-25 prikazan je primer periodičnog ponašanja, u vidu odgovarajuće vremenske serije (Slika 7-25a) i faznog portreta (Slika 7-25b), za vc=1. Usamljeni pik u Furijeovom spektru snage potvđuje oscilatorno ponašanje (slika 7-26a), dok maksimalni Ljapunovljev eksponent konvergira ka X=-0,0296, što takođe ukazuje na periodično ponašanje (slika 7-26b). U svim izvedenim simulacijama, tranzijentno ponašanje je zanemareno.\r\nSlika 7-25. Vremenska serija U(t)(a) i odgovorajući fazni portret (b) za vc=1 (oscilacije).\r\nSlika 7-26. Furijeov spektar snage (a) i određivanje maksimalnog Ljapunovljevog\r\neksponenta (b) za vremensku seriju U (t) prikazanu na slici 7-25a: (a) usamljeni pik u spektru snage ukazuje na oscilatorno ponašanje modela; (b) najveći Ljapunovljev eksponent konvergira ka X= - 0,0296.\r\nSada prelazimo na analizu sistema (7.6), za različite vrednosti parametra jačine trenja c, i vremenskog kašnjenja T. Uvodeći samo parametar c, pri čemu je T=0, i variranjem parametra vc, ne dolazi do novog tipa dinamičkog ponašanja modela, izuzev za jako niske vrednosti parametra vc (<10-8), kada blok prelazi u ravnotežno stanje.\r\nSada uvodimo vremensko kašnjenje T u izrazu za trenje, i posmatramo dinamičko ponašanje modela. Za konstantno c=1, ne dolazi do promena u ponašanju modela, što je i očekivano, s obzirom na to da nema promene u jačini trenja [Burridge, Knopoff, 1967; Carlson, Langer, 1989]. Međutim, ukoliko menjamo vrednost parametra c, onda uvedeno vremensko kašnjenje uzrokuje pojavu ,,oscilatorne smrti'' [Yamaguchi, Shimizu, 1984; Aronson i dr., 1990, Reddy i dr., 1998]. Ova pojava je vrlo značajna i interesantna sa inženjerske tačke gledišta, jer ukazuje na činjenicu da sa povećanjem parametra T može doći i do smirivanja kretanja bloka (a time i kretanja duž raseda), što može dovesti do aseizmičkog kretanja duž raseda bez pojave potresa.\r\nSa povećanjem vrednosti parametra T, blok prelazi iz ravnotežnog stanja u osiclatorno kretanje, i opet u ravnotežno stanje. Na slici 7-27 prikazan parametarski prostor (T,Vc) sistema sa jednim blokom, na kojem su izdvojene oblasti različitog dinamičkog ponašanja (ravnotežno stanje i oscilatorno kretanje), za fiksiranu vrednost parametra c (c=3).\r\nSlika 7-27. Parametarski prostor (T,Vc) sistema sa jednim blokom, za c=3. ES i PM predstavljaju redom skraćenice za ravnotežno stanje i oscilacije bloka. Dijagram je konstruisan za veličinu koraka 0,1, i za T i za vc. Vremenske serije i fazni portreti za tačku 1 i 2 prikazani su na Slici 7-28.\r\nPrimeri ravnotežnog stanja i oscilacija za tačke 1 i 2 na slici 7-27, redom, prikazani su u obliku odgovarajućih vremenskih serija i faznih portreta na slici 7-28.\r\nSlika 7-28. Vremenska serija U (t) i odgovarajući fazni portret za: (a) ravnotežno stanje (tačka 1 sa slike 7-27; vc=7, T=3, c=3), (b) periodično ponašanje (tačka 2 sa slike 7-27; vc=7, T=4, c=3). Fazni portret na slici (b) dobijen je zanemarujući tranzijentno ponašanje.\r\nTakođe, dve odvojene oblasti sa različitim dinamičkim ponašanjem (ravnotežno stanje i oscilacije) pojavljajuju se i kada se vremensko kašnjenje T drži konstantnim, a parametar c se menja. Na slici 7-29 prikazan je ovaj slučaj, kada je T=3.\r\nSlika 7-29. Parametarski prostor (c,vc) sistema sa jednim blokom, za T=3. Dijagram je konstruisan za veličinu koraka 0,1, i za c i za vc. ES i PM predstavljaju redom skraćenice za ravnotežno stanje i oscilacije bloka.\r\nSada osmatramo ponašanje modela istovremeno varirajući oba uvedena parametra, c i T, za vc=1 [Vieira, 1999]. U ovom slučaju, takođe se izdvajaju dve oblasti sa različitim dinamičkim ponašanjem - ravnotežno stanje i oscilacije (slika 7-30). Analiza pokazuje da se periodično ponašanje osmatra za vrednosti parametra c < 10.\r\nSlika 7-30. Parametarski prostor (c,T) sistema sa jednim blokom, za vc=1. Dijagram je konstruisan za veličinu koraka 0,1, i za T i za c. ES i PM predstavljaju redom skraćenice za ravnotežno stanje i oscilacije bloka.\r\n7.2.1.2. MODEL SA DVA BLOKA\r\nSlična analiza izvedena je i za model sa dva bloka (7.7). U svim ispitivanim slučajevima, usvaja se da je k=k2=k= 1. U prvom koraku, vrši se analiza samo početnog modela (c=c2=i) bez vremenskog kašnjenja (T=0), i sa fiksiranim vrednostima parametara vci=vc2=1. U ovom slučaju, ispitivani model pokazuje samo oscilatorno ponašanje. Odgovarajuća vremenska serija i fazni portret za oscilatorno ponašanje prikazani su na slici 7-31. Usamljeni pik u Furijeovom spektru snage potvrđuje oscilatorno ponašanje modela (slika 7-32a), dok negativna vrednost maksimalnog Ljapunovljevog eksponenta X= -0,039 takođe ukazuje na periodično ponašanje (Slika 7- 32b).\r\nSlika 7-31. Vremenska serija U(t) (a) i odgovarajući fazni portret (b) za c=c2=1, T=0, k=k2=k=1 i vci=vc2=1 (oscilacije).\r\nU sledećem koraku, najpre se uvode parametri jačine trenja c i c2. Istovremenim variranjem ova dva uvedena parametra, sa konstanom vrednošću parametara vci=vc2=1, model pokazuje kvaziperiodično ponašanje. Jedino za jednake vrednosti parametara c i c2, model pokazuje oscilatorno ponašanje.\r\nSlika 7-32. Furijeov spektar snage (a) i određivanje maksimalnog Ljapunovljevog eksponenta (b) za vremensku seriju U (t) prikazanu na slici 7-31a: (a) usamljeni pik u spektru snage ukazuje na oscilatorno ponašanje modela; (b) maksimalni Ljapunovljev eksponent konvergira ka negativnoj vrednosti X= - 0,039.\r\nSlično, ukoliko se uvede vremensko kašnjenje u izrazu za trenje za oba bloka, tj i t2, sa identičnom jačinom trenja ci=c2, model takođe pokazuje kvaziperiodično ponašanje, izuzev u slučaju jednakih vrednosti parametara t] = T2, kada je ponašanje periodično. Primer kvaziperiodičnog ponašanja, za vrednosti parametara T=4, t2=5, ci=c2=2, i vci=vc2=\\ prikazan je na slici 7-33, u obliku odgovarajuće vremenske serije i faznog portreta. Dva pika u Furijeovom spektru snage potvrđuju kvaziperiodično ponašanje modela (slika 7-34a), dok je druga potvrda kvaziperiodičnosti data preko približno nulte vrednosti maksimalnog Ljapunovljevog eksponenta (slika 7-34b).\r\nSlika 7-33. Vremenska serija U^(t) (a) i odgovarajući fazni portret (b) za TJ=4, T2=5, cj=c2=2, k]=k2=k=1 i vcj=vc2=1 (kvaziperiodično kretanje).\r\nSlika 7-34. Furijeov spektar snage (a) i određivanje maksimalnog Ljapunovljevog\r\neksponenta (b) za vremensku seriju U(t) prikazanu na slici 7-33a: (a) dva pika u spektru snage ukazuju na kvaziperiodično ponašanje; (b) maksimalni Ljapunovljev eksponent približno konvergira ka A~0.\r\nU sledećem koraku, vrši se razmatranje kretanja blokova za slučaj kada se prvi blok drži u ravnotežnom stanju, za c^=0,1 i t^=3,0 dok se parametri c2 i T2 za drugi blok variraju. Rezultati analize su prikazani na slici 7-35, gde se jasno može uočiti da preovladava kvaziperiodično ponašanje, sa dve relativno male oblasti ravnotežnog stanja (u parametarskom prostoru) za male vrednosti parametra jačine trenja (c2<1). Analiza takođe pokazuje da kvaziperiodično ponašanje preovladava za c2< 10, bez pojave drugih tipova ponašanja. Oscilacije blokova (periodično ponašanje) javljaju se samo kao tranzijentna pojava, sa vrlo brzim prelaskom u kvaziperiodično ponašanje.\r\nSlika 7-35. Parametarski prostor (c2,T2) sistema sa dva bloka, za konstantne vrednosti cj=0,1, 17=3,0 (prvi blok u ravnotežnom stanju). Dijagram je konstruisan za veličinu koraka 0,1, i za T2 i za c2. ES i QM predstavljaju redom skraćenice za ravnotežno stanje i kvaziperiodično ponašanje. Vrednosti ostalih parametara su kao na slici 7-33.\r\nMeđutim, ukoliko prvi blok pokazuje oscilatorno ponašanje (slika 7-30), a parametri c2 i T2 za drugi blok se menjaju, posmatrani model sa dva bloka pokazuje uglavnom kvaziperiodično ponašanje, pri čemu se periodično ponašanje i ravnotežno stanje javljaju samo za izolovane vrednosti u parametarskom prostoru sistema sa dva bloka (slika 7-36).\r\nSlika 7-36. Parametarski prostor (c2,T2) sistema sa dva bloka, za konstantne vrednosti parametara za prvi blok c7=0,2 ,t7=0,5 (oscilatorno ponašanje za prvi blok). Dijagram je konstruisan za veličinu koraka 0,1, i za T2 i za c2. Crne tačke označavaju oscilatorno ponašanje, sive tačke - ravnotežno stanje, belo područje - kvaziperiodično ponašanje. Vrednosti ostalih parametara su: vc7=vc2=1, k7=k2=k=1.\r\nPored periodičnog i kvaziperiodičnog ponašanja, utvrđena je i pojava ponašanja nalik determinističkom haosu, ali samo kao tranzijentnog fenomena (slika 7-37). Ovde se termin ,,tranzijentni'' koristi da bi se naglasio relativno kratak vremenski interval u kome se ovo ponašanje opaža u odnosu na celokupno trajanje analizirane vremenske serije (mala slika iznad slike 7-37a). Kontinualni širokopojasni šum u spektru snage (slika 7-38) potvrđuje ponašanje posmatranog modela nalik haotičnom.\r\nSlika 7-37. Vremenska serija Ui (t) (a) i odgovarajući fazni portret (b) za ci=0.2, T^=0.5 (oscilatorno ponašanje za prvi blok), c2=0.2 i T2=3 (deterministički haos). Vrednosti ostalih parametara su: vf =v2c = 1, k^=k2=k= 1. Celokupna vremenska serija, za 2000 jedinica, prikazana je na maloj slici iznad slike 7-37a, na kojoj je ,,haotična'' oblast izdvojena pravougaonikom.\r\nSlika 7-38. Širokopojasni šum u Furijeovom spektru snage može ukazivati na postojanje ponašanja nalik relativno slabom tranzijentnom determinističkom haosu.\r\nNeophodno je naglasiti da je određivanje vrednosti najvećeg Ljapunovljevog eksponenta, za slučaj tranzijentnog ponašanja nalik haotičnom, još uvek u fazi istraživanja. U opštem slučaju, ,,proglašavanje'' ponašanja tranzijentnim haosom moguće je samo na sonovu vrednost Ljapunovljevog eksponenta u konačnom vremenskom intervalu (Stefanski, i dr. 2010). U okviru ovog istraživanja, pojava ponašanja nalik determinističkom haosu je potvrđena pozitivnom vrednošću maksimalnog Ljapunovljevog eksponenta, koristeći dve metode (Wolf, i dr., 1985; Rosenstein i dr., 1993). U oba slučaja, vrednosti maksimalnog Ljapunovljevog eksponenta su pozitvne i istog reda veličine: AMAX=0,0016 i 0,0019, redom (slika 7-39).\r\nSlika 7-39. Određivanje vrednosti maksimalnog Ljapunovljevog eksponenta za vremensku seriju prikazanu na slici 7-37 (2): (a) AMAX=0,0016 (Wolf et al., 1985); (b) Amax~0,0019 (Rosenstein et al., 1993). Nagib pravih linija koji označava rast efektivne brzine razdvajanja S(An) predstavlja robusnu ocenu vrednosti maksimalnog Ljapunovljevog eksponenta. Rezultati su dobijeni za 1000 referentnih tačaka i malu okolinu e=0,1-0,15.\r\nMeđutim, kako sistem (7.7) predstavlja sistem diferencijalnih jednačina sa kašnjenjem, to znači da su početni uslovi dati u obliku početne funkcije na intervalu [-T,0]. Prema tome, pristupilo se izvođenju dodatne analize pojave ponašanja nalik haotičnom usrednjavanjem po početnim uslovima, u intervalu U\" =[0;0,003], U20 = [0;0,05],\r\nUj° = [0;0,003], U20 = [0;0,07]. U svim ispitivanim slučajevima maksimalni Ljapunovljevi eksponenti kovergiraju ka pozitivnim vrednostima za različite početne uslove (slika 7-40). Štaviše, ove vrednosti su istog reda veličine (10-3) kao i kod prvobitno ispitivanog slučaja (slika 7-39).\r\nSlika 7-40. Određivanje vrednosti maksimalnog Ljapunovljevog eksponenta, po postupku iz (Wolf, i dr., 1985), usrednjavanjem po početnim uslovima u intervalu:\r\nU° =[0,0.003], U20 =[0;0,05], U° =[0;0,003], U20 =[0;0,07]. U svim ispitivanim slučajevima, maksimalni Ljapunovljev eksponent konvergira ka pozitivnim vrednostima reda veličine 10-3 (kao i za slučaj prikazan na slici 7-39).\r\nU poslednjem koraku, pristupilo se analizi uticaja različitih početnih uslova na dinamičko ponašanje modela sa dva bloka. Izvedena analiza ukazuje na to da je dinamičko ponašanje ovog modela zavisno od početnih uslova, odnosno, za različite početne uslove, u okolini ravnotežnog stanja (0;0,1), u rasponu vrednosti [-1,1], dolazi do pojave različitog dinamičkog ponašanja, odnosno kretanje dva bloka prelazi iz ravnotežnog stanja u oscilacije. Slika 7-41 pokazuje atraktore kada se vrednosti početnih uslova menjaju oko ravnotežne tačke sistema (U° = U20 = 0, U° = U20 = 0,1).\r\nKao što se vidi sa slike 7-41, sa većim udaljenjem od ravnotežnog stanja, sistem pokazuje oscilatorno ponašanje, ukazujući na moguću pojavu bistabilnosti u sistemu. Ova osetljivost modela bloka sa oprugom na promenu početnih uslova već je osmatrana u radu [Szkutnik i dr., 2003]. Njihova analiza je pokazala da karakter kretanja u modelu sa tri bloka, za određenu vrednost kontrolnih parametara, zavisi od početnih uslova. Drugim rečima, samo promenom početnog položaja jednog od blokova, dolazi do prelaza iz kvaziperiodičnog i nesinhronizovanog kretanja u periodično kretanje sa sinhronizacijom dva krajnja bloka.\r\nSlika 7-41. Prostor početnih uslova (U° = U20 = U0, Uj0=U20=U°) sa pojavom ravnotežnog stanja (ES) i periodičnog ponašanja (PM), za vrednosti parametara:\r\nc7=c2=0,1, TI= T2=3,0. Dijagram je konstruisan za veličinu koraka 0,1 za U0 i U, redom. Vrednosti ostalih parametara su: vf = v2c = 1, k7=k2=k=1.\r\n7.2.2. ANALIZA PONAŠANJA REŠENJA BEKEROVOG SISTEMA ZA MODEL SA JEDNIM BLOKOM\r\nZa razliku od Karlson-Langer sistema, analiziranog u prethodnom odeljku, Beker [2000] polazi od Diterihovog konstitutivnog zakona trenja sa dve promenljive stanja, 6, i 62:\r\ngde su parametri A, B^ i L konstante; V* i T* označavaju referentnu brzinu i čvrstoću na smicanje, redom, dok 6^ i V zavise implicitno od vremena t. Parametar LT odgovara kritičnom rastojanju DC koje blok pređe za vreme kretanja. Promena varijable stanja 6 definiše se na sledeći način:\r\nšto je prvobitno predloženo od strane Ruine [1983].\r\nPromena trenja u zavisnosti od početne brzine Vo (pod dejstvom sile), brzine kretanja bloka V i krutosti opruge, K>0, data je u sledećem obliku:\r\nčime je kompletiran matematički model bloka sa oprugom (Slika 7-42). Primetimo da u ovom slučaju blok nije povezan preko opruge sa pokretnom pločom, već je zadnji kraj opruge izložen dejstvu sile.\r\nSlika 7-42. Blok sa oprugom sa jednim stepenom slobode. Blok se pokreće konstantnom brzinom V0, dok se trenje duž kontaktne površi menja po zakonu zavisnom od brzine i dve promenljive stanja.\r\nJednačine (7.13), (7.14) i (7.15) mogu biti predstavljene i u bezdimenzionalnom obliku, uvođenjem sledećih promenljivih:\r\nkada se dobija rezultujući sistem jednačina u bezdimenzionalnom obliku:\r\npri čemu je početna brzina V° izjednačena sa brzinom V* bez gubitka opštosti. Sistem (7.17) opisuje kretanje bloka u funkciji bezdimenzionalne brzine x, napona (stresa) y, i promenljive stanja z.\r\nPrethodna istraživanja sistema jednačina (7.17) pokazala su da za fiksirane vrednosti parametara fii=1, ^2=0,84 i p=0,048, i pri promeni samo jednog kontrolnog parametra K (bezdimenzionalna krutost opruge), dolazi do pojave determinističkog haotičnog ponašanja za K<0.85 [Gu, i dr., 1984; Becker, 2000].\r\nU ovom slučaju, pretpostavlja se periodična perturbacija parametara K i fi^ u obliku pozitivnih periodičnih funkcija vremena:\r\ngde Sp, SK, rnp, i rnK predstavljaju konstantne amplitude oscilacija i ugaone frekvencije, redom. Pri promeni vrednosti parametra, usvojeno da je Sp0,85 (u suprotnosti sa rezultatima Becker, 2000) što je potvrđeno širokopojasnim,\r\nkontinualnim šumom u Furijeovom spektru snage (slika 7-44a) i pozitivnom vrednošću maksimalnog Ljapunovljevog eksponenta (slika 7-44b).\r\nSlika 7-43. Vremenska serija x(t) (a) i odgovarajući fazni portret (b) sistema (7.14) za za amplitude perturbacije Sp=0,9 i SK=0,5, dok su vrednosti ostalih parametara (za sistem koji u neperturbovanom stanju ispoljava periodične oscilacije) održavane konstantnim: fii=1, k=1,5, ^2=0,84, r=0,048, ^g=0,8, wK=0,5 (,,slabo'' determinističko haotično ponašanje).\r\nSlika 7-44. Furijeov spektar snage (a) i određivanje maksimalnog Ljapunovljevog eksponenta (b) za vremensku seriju x(t) prikazanu na slici 7-43a: (a) kontinualni širokopojasni šum u spektru snage ukazuje na slabo determinističko haotično ponašanje; (b) maksimalni Ljapunovljev eksponent približno konvergira ka A~0,04.\r\n7.2.3. ANALIZA DINAMIKE MADARIAGA SISTEMA ZA MODEL SA JEDNIM BLOKOM\r\nDrugi model predstavljen je sistemom jednačina predloženim od strane Madariage [Erickson i dr., 2008]. Ove jednačine kretanja bloka, sa Diterih-Ruina zakonom trenja, u originalnom obliku date su na sledeći način:\r\ngde je M masa bloka, a k konstanta opruge, koja odgovara elastičnim svojstvima stenske mase u neposrednoj blizini rasedne zone [Scholz, 2002]. Prema [Dieterich, Kilgore, 1994] parametar L odgovara rastojanju koje blok pređe za vreme kretanja. Parametri A i B predstavljaju empirijske konstante, koje zavise od svojstava stenske mase. Prema [Rice, 1993: Rice i dr., 2001], parametar A meri direktnu zavisnost trenja od brzine (tzv.,,direktni efekat\"), dok (A-B) predstavlja meru zavisnosti trenja od brzine u stacionarnom stanju (kada se blok kreće). Radi lakše analize, pogodno je izvršiti bezdimenzionalizaciju sistema (7.19) uvodeći nove promenljive 6', v', u' i t' na sledeći način: 6=A6', v=v°v', u=Lu', t=(L/v°)t', a potom vraćanjem na prethodne promenljive 6, v, u i t, čime dobijamo sistem jednačina u sledećem obliku:\r\ngde s=(B A)/A predstavlja osetljivost bloka na promenu brzine gornje ploče, £ =(kL)/A\r\n1/2\r\nje bezdimenzionalna konstanta opruge, a y= (k/M) (L/v°) predstavlja bezdimenzionalnu frekvenciju trzanja [Erickson i dr., 2008].\r\n7.2.3.1. PRIMENJENE METODE\r\nS obzirom na logaritamski član u sistemu jednačina (7.20), analitičko rešavanje ovog sistema je jako otežano. Drugim rečima, logaritamski član čini posmatrani sistem numerički jako ,,krutim'' za uobičajene vrednosti parametara, u smislu da je potreban jako mali korak iteracije kako bi numerička integracija bila uspešno izvedena. Ovaj problem se može prevazići pojednostavljenjem posmatranog sistema, na primer regularizacijom nelinearnog člana za trenje za brzine bliske nuli [Lapusta, Rice, 2003; Szkutnik i dr., 2003]. Međutim, vrlo često ovakve aproksimacije su uzrokovale nemogućnost osmatranja kompleksnijeg ponašanja u posmatranim sistemima kada se parametrima dodele realne vrednosti. Iz tog razloga, potrebno je odabrati numeričku shemu koja će dozvoliti primenu punog nelinearnog člana. Shodno tome, umesto primene neke od eksplicitnih metoda koje mogu postati nestabilne u jednom trenutku integracije, ovde se usvaja metoda zadnjih razlika (ili metoda podeljenih razlika unazad) (BDF)4, koja pripada klasi implicitnih postupaka integracije, a koja je već primenjena pri rešavanju numerički krutih sistema diferencijalnih jednačina [Erickson i dr., 2008]. Ova metoda pripada grupi linearnih višekoračnih metoda numeričkog rešavanja običnih diferencijalnih jednačina [Mayers, 2003].\r\nŠto se tiče analize lokalnih bifurkacija, bifurkacione krive su prvo dobijene analitički, a potom su rezultati upoređeni sa rezultatima dobijenim numeričkim postupkom, primenom softverskog paketa DDE-BIFTOOL, koji predstavlja skup MATLAB alatki za numeričku bifurkacionu analizu sistema diferencijalnih jednačina sa kašnjenjem [Engelborghs, 2000; Engelborghs i dr., 2001]. DDE-BIFTOOL je do sada uspešno primenjivan u oblasti biologije, hemije i fizike [Luyzanina, i dr., 2001; Haegeman i dr., 2002]. Numerička proračunavanja (simulacije) su izvedena na osnovu algoritma u programskom jeziku C++ (mr Igor Franović, dipl.fiz.).\r\nU cilju određivanja vrednosti maksimalnog Ljapunovljevog eksponenta, korišćen je softver M. Perca sa Univerziteta u Mariboru, koji je razvijen za potrebe analize\r\nnelinearnih vremenskih serija [Perc, 2012]. Ovaj softver se oslanja na dobro poznati Vulfov metod, opisan u poglavlju ,,Metodologija istraživanja''.\r\n7.2.3.2. ANALIZA PONAŠANJA REŠENJA U USLOVIMA PERIODIČNE PERTURBACIJE KONTROLNIH PARAMETARA\r\nU prvoj fazi, pristupa se razmatranju dinamike ponašanja modela pri periodičnoj promeni parametara e and £\r\ngde su e(t) i £(t) pozitivne periodične funkcije vremena:\r\ni SE0,34 i S^>0,36. Naglasimo da za S^>0,45, sistem (7.21) postaje numerički krut, u smislu da su potrebni jako mali koraci integracije da bi se došlo do rešenja (reda\r\nveličine 10-5).\r\nU sledećem koraku, pristupamo analizi dinamike sistema (7.21) pri varijaciji oba parametra £ i £ Kao i u prethodnom slučaju, za vrednosti ugaonih frekvencija izabrane su one koje odgovaraju oscilacijama bloka u neprturbovanom stanju, ali haotičnom ponašanju pri periodičnoj promeni samo jednog parametra. Dinamika sistema se sada osmatra samo za promenu amplituda SE i S^ (Slika 7-58).\r\nSlika 7-58. Parametarski prostor (SE,S^) sistema (7.21), pri frekvenciji oscilacija koje su bliske frekvenciji oscilacija modela u neperturbovanom stanju, ali za haotično ponašanje pri periodičnoj varijaciji samo jednog od parametara (w£=0,2; we=0,3). Dijagram je konstruisan za veličinu koraka 0,05, i za SE i za S^. Za svaki korak, parametri sistema uzimaju konstantne vrednosti koje odgovaraju stanju sistema blizu ravnotežne tačke, ali za pojavu grničnog ciklusa: e=0.4, ^=0.5 i y=0.80. QP i C predstavljaju skraćenice za kvaziperiodično ponašanje i deterministički haos, redom. Periodično ponašanje se osmatra samo za jako niske vrednosti amplituda (<10-3). Odgovarajuće vremenske serije i fazni portreti za tačke 1 i 2 dati su na slici 7-59.\r\nKao što je jasno sa slike 7-58, u ovom slučaju do pojave determinističkog haotičnog ponašanja dolazi za mnogo niže vrednosti amplituda perturbacija, SE i S^, u poređenju sa amplitudama perturbacija samo jednog parametra e ili £ Tipične vremenske serije i fazni portreti za kvaziperiodično ponašanje (tačka 1 na slici 7-58) i haotično ponašanje (tačka 2 na slici 7-58) prikazani su na slici 7-59. Dva izražena pika u spektru snage (slika 7-60a) potvrđuju kvaziperiodično ponašanje, dok kontinualni širokopojasni šum u spektru snage (slika 7-60b) potvrđuje pojavu determinističkog haosa.\r\nSlika 7-59. Vremenske serije V(t) i odgovarajući fazni portreti za tačke 1 i 2 na slici 7- 58, za sledeće vrednosti parametara: (1) 4=0,1; 4=0,05 (kvaziperiodično ponašanje); (2) 4=0,2; 4=0,2 (deterministički haos). U svakom koraku, parametri sistema se drže konstantnim za vrednosti blizu fiksne (singularne) tačke, ali za koje blok ispoljava oscilatorno ponašanje, e= 0,4; ^=0,5 i 7=0,8. Frekvencije oscilacija su bliske frekvenciji oscilacija modela u neperturbovanom stanju, ali za haotično ponašanje pri periodičnoj varijaciji samo jednog od parametara (w£=0,2; w£=0,3).\r\nSlika 7-60. Furijeov spektar snage za vremenske serije prikazane na slici 7-59: (a) dva izražena pika u spektru snage za vremensku seriju na slici 7-59(1) ukazuju na pojavu kvaziperiodičnog ponašanja u modelu; (b) širokopojasni kontinualni šum u Furijeovom spektru snage za vremensku seriju sa slike 7-59(2) ukazuje na determinističko haotično ponašanje u modelu.\r\nDeterministički haos je dalje potvrđen pozitivnom vrednošću maksimalnog Ljapunovljevog eksponenta, koristeći dva različita postupka (Wolf, i dr., 1985; Rosenstein, i dr., 1993), kao što je prikazano na slici 7-61.\r\nSlika 7-61. Određivanje vrednosti maksimalnog Ljapunovljevog eksponenta za vremensku seriju prikazanu na slici 7-59(2): (a) LMAX=0,009 (Wolf et al., 1985); (b) XMAX~0,009 (Rosenstein et al., 1993). Nagib pravih linija koji označava brzinu ekspanzije S(An) prvobitno susednih tačaka predstavlja robusnu ocenu vrednosti maksimalnog Ljapunovljevog eksponenta. Rezultati su dobijeni za 1000 referentnih tačaka i malu okolinu £=0,1-0,5.\r\n7.2.3.3.STANDARDNA LOKALNA BIFURKACIONA ANALIZA MODELA SA VREMENSKIM KAŠNJENJEM\r\nU drugom delu analize, u izraz za trenje u sistemu jednačina (7.21) uvodi se vremensko kašnjenje kojim se simulira odloženi efekat rudarskih radova na generisanje potresa duž postojećih raseda, čime se dobija sledeći sistem diferencijalnih jednačina:\r\nSistem (7.23) ima samo jedno stacionarno rešenje (6,u,v)=(0,0,1), koje odgovara stacionarnom kretanju bloka. Linearizacijom sistema (7.23) i zamenom 0=AeAt, u=BeAt, v=CeAt i v(t-r)=CeA(t~x) dobijamo sistem algebarskih jednačina po A, B i C. Sistem ima netrivijalno rešenje ukoliko je Jakobijan sistema (7.23) jednak nuli:\r\nodnosno, ukoliko je karakteristična jednačina sistema (7.23):\r\n7.2.3.3.1. STABILNOST SISTEMA ZA SVE VREDNOSTI VREMENSKOG KAŠNJENJA\r\nPre početka analize, neophodno je odrediti deo parametarskog prostora za koji je stacionarno rešenje stabilno za sve vrednosti T. U tom cilju, pogodno je primeniti Rušeovu teoremu [Titchmarsh, 1939], koja je opisana u poglavlju ,,Metodologija istraživanja''. Prema ovoj teoremi, neophodno je odrediti vrednosti kontrolnih parametara tako da na konturi C važi:\r\ngde su\r\nUkoliko pođemo od sistema (7.24) dobijamo:\r\nStoga, odavde sledi:\r\nS druge strane, za y(k), dobijamo:\r\nodnosno:\r\npretpostavljajući da je >na konturi C. Iz prethodnog sledi:\r\nPrema tome, uslovi Rušeove teoreme su ispunjeni ukoliko je:\r\nStoga, ukoliko je uslov (7.28) zadovoljen, onda je stacionarno rešenje stabilno za bilo koju vrednost vremenskog kašnjenja T. Za one delove parametarskog prostora, za koje ne važi odnos (7.28), stabilnost stacionarnih rešenja za vremensko kašenjenje različito od nule mora se dalje istraživati, određivanjem bifurkacionih krivih koje ,,razdvajaju'' delove parametarskog prostora sistema sa različitim dinamičkim ponašanjem.\r\n7.2.3.3.2. LOKALNA STABILNOST I BIFURKACIJE STACIONARNOG REŠENJA\r\nJednačina (7.25) može se zapisati u sledećem obliku:\r\nKao što je već opisano u poglavlju ,,Metodologija istraživanja'' od interesa za stabilnost sistema (7.23) jesu sopstvene vrednosti (gde je i=1,2,3...) karakteristične jednačine (7.29). Ukoliko je X=0, iz (7.26) sledi da je rešenje +/- y=0. Međutim, ukoliko je X=irn, iz (7.29) dobijamo:\r\nIz (7.30) možemo izdvojiti dve posebne jednačine za realni (A) i imaginarni deo (B):\r\nKvadriranjem prve i druge jednačine sistema (7.31) i njihovim sabiranjem, a na osnovu Ojlerove teoreme, dobija se:\r\nodnosno:\r\nšto predstavlja izraz za parametar s u funkciji ostalih kontrolnih parametara sistema, £ i 7.\r\nS druge strane, iz (7.30) dalje sledi da je:\r\nodakle sledi izraz za parametar £ u funkciji ostalih kontrolnih parametara sistema, e, y i T:\r\nI, konačno, deljenjem druge jednačine sa prvom u sistemu (7.31) dobijamo da je:\r\nodnosno dobija se izraz za T u funkciji £ i 7:\r\ngde je c bilo koji nenegativni ceo broj takav da je TC > 0.\r\nIako sam oblik rešenja karakteristične jednačine ukazuje na moguću pojavu Hopfovih bifurkacija, strogi dokaz ove tvrdnje zahteva duži analitički postupak, i neće biti ovde razmatran [Belair, Campbell, 1994; Wiggins, 2000; Kuznetsov, 2004].\r\nDovoljno je reći da se prethodno date parametarske jednačine za e, £ i T poklapaju sa Hopfovim bifurkacionim krivama, prikazanim na slikama 7-62 i 7-63. Na slici 7-62 prikazane su bifurkacione krive T(S) za fiksirane vrednosti parametara £=0,5 i 7=0,8. Na dijagramu, datom na slici 7-63, prikazane su bifurkacione krive £(e), za fiksirane vrednosti parametara T=12 i 7=0.8.\r\nSlika 7-62. Hopfove bifurkacione krive T(Č), za fiksirane vrednosti parametara £ = 0,5, i 7=0,8. Predznak +/- redom označava destabilizirajuću (direktnu) ili stabilizirajuću (inverznu) Hopfovu bifurkaciju, sa povećanjem vrednosti vremenskog kašnjenja T. Odgovarajuće vremenske serije i fazni portreti za tačke a, b, c i d prikazani su na slici 7- 64.\r\nU sistemu (7.23) dolazi do pojave i superkritičnih i subkritičnih Hopfovih bifurkacija, što nije od posebnog značaja za ovu analizu. Mnogo je značajnije to da li su ove bifurkacije direktnog ili inverznog tipa [Belair, Campbell, 1994], što rezultuje u stvaranju/poništavanju nestabilne ravni u prostoru stanja sistema onog trenutka kada se odgovarajuća kriva pređe u određenom pravcu. Drugim rečima, direktna (inverzna) bifurkacija ima destabilizirajući (stabilizirajući) uticaj na dinamiku sistema.\r\nSlika 7-63. Hopfove bifurkacione krive £(e) za fiksirane vrednosti parametara T=12 i 7=0,8.\r\nS obzirom na to da je glavna pažnja analize usmerena na uticaj vremenskog kašnjenja T na dinamiku sistema, pristupilo se bližem ispitivanju bifurkacija sa slike 7- 62 koje se javljaju u sistemu sa povećanjem T. U tom cilju, potrebno je bliže odrediti promenu realnog dela korena karakteristične jednačine (7.25), za kritične vrednosti parametra T:\r\nDrugim rečima, ukoliko je dRe(A)/dr > 0, onda se broj korena sa pozitivnim realnim delom uvećava (destabilizirajuća bifurkacija), dok se njihov broj smanjuje ukoliko je prethodno pomenuti izvod negativan (stabilizirajuća bifurkacija). Bifurkacione krive na slici 7-62 prikazane su sa predznakom +/- kako bi bio naglašen stabilizirajući/destabilizirajući uticaj povećanja/smanjenja vremenskog kašnjenja. U ovom slučaju, diferenciranje se vrši na sledeći način. Polazeći od izvoda karakteristične jednačine A(A) = 0, korišćenjem pravila izvoda složene funkcije, dobijamo\r\nšto se može predstaviti i u sledećem obliku:\r\npod uslovom da važi\r\nTada, smenom l=irn i izdvajanjem realnog člana, konačni izraz za željeni izvod je oblika\r\nUkoliko analizu bifurkacionih krivih počnemo od prve koja se javlja sa povećanjem vremenskog kašnjenja, dolazimo do zaključka da su odgovarajući izvodi duž bifurkacione krive negativni, što ukazuje na to da sistem prolazi kroz inverznu bifurkaciju. Drugim rečima, stanje sistema se menja od oscilatornog ka stabilnom stacionarnom stanju. Radi poređenja, na slici 7-64 prikazane su vremenske serije i fazni portreti promenljive v za fiksirane vrednosti parametara e, £ y i T. Ova pojava, kada ravnotežno stanje ,,vraća'' svoju stabilnost sa povećanjem T predstavlja primer ,,oscilatorne smrti'' izazvane vremenskim kašnjenjem [Campbell, 2007], što predstavlja predmet istraživanja u različitim oblastima nauke [Prasad i dr., 2010; Reddy i dr., 1998]. Pojava oscilatorne smrti naglašava potencijalno zanimljivu ulogu malih vrednosti vremenskog kašnjenja, koje može biti uzrok smanjenja seizmičke aktivnosti. S druge strane, sa slike 7-63 jasno je da dalje povećanje T ima destabilizirajući efekat na rasednu dinamiku. Ne samo da vremensko kašnjenje uzrokuje promenu stabilnosti stacionarnog stanja, izazivajući na taj način pojavu periodičnih oscilacija, već može doći do pojave niza bifurkacija, koje, konačno, vode ka haotičnom ponašanju. Preciznije, stanje sistema se menja prateći kvaziperiodični put u haos (Ruel-Takens-Njuhaus put u haos5) [Ruelle, Takens, 1971; Newhouse i dr., 1978], tako da nakon dve superkritične Hopfove bifurkacije, stacionarno kretanje postaje vrlo nestabilno, uzrokujući dinamiku sistema, koja, u konačnom vremenskom intervalu, konvergira ka stranom atraktoru. Dakle, samo povećavajući vrednosti parametra T, tj. za vrednosti parametra T: 0, 10, 13 i 20, redom, i samo blagom promenom vrednosti ostalih parametara, dinamika bloka se menja od ravnotežnog stanja (fiksne tačke), preko periodičnih oscilacija - graničnog ciklusa (prva Hopfova bifurkacija) i kvaziperiodičnih oscilacija - torusa (druga Hopfova bifurkacija) do determinističkog haosa.\r\nKao što je već jasno iz prethodne analize, u ovom slučaju takođe, kao i pri periodičnoj perturbaciji parametara, do pojave determinističkog haosa dolazi za značajno manje vrednosti parametra s (0,5) u poređenju sa vrednostima u [Erickson i dr., 2008].\r\nPut ka determinističkom haosu dalje je potvrđen izračunavanjem Furijeovog spektra snage za oscilacije, torus i haotične orbite, koji su prikazani na slici 7-65 i 7-66, za odgovarajuće vremenske serije na slikama 7-64b, 7-64c i 7-64d, redom. Usamljeni pik u spektru snage na slici 7-65a ukazuje na oscilatorno ponašanje sistema, dok dva pika na slici 7-65b ukazuju na kvaziperiodično ponašanje na torusu. Kontinualni širokopojasni šum u spektru snage na slici 7-66 ukazuje na prisustvo stranog atraktora.\r\nSlika 7-64. Vremenska serija promenljive v i odgovarajući fazni portreti za vrednosti parametara (a) T=0, e = 0,2, £ = 0,5 i y = 0,8 (fiksna tačka); (b) T=10, e = 0,3, £ = 0,5 i Y=0,8 (oscilacije); (c) T=13, e = 0,5, £ = 0,5 i y = 0,8 (torus); (d) T=20, e = 0,5, £ = 0,5 i Y= 0,8 (deterministički haos).\r\nSlika 7-65. (a) Usamljeni pik u spektru snage ukazuje na oscilatorno ponašanje modela. (b) Dva pika u spektru snage ukazuju na pojavu kvaziperiodičnog ponašanja na torusu (druga Hopfova bifurkacija). Vrednosti parametara su iste kao na slikama 7-64b i 7-64c, redom.\r\nSlika 7-66. Kontinualni širokopojasni šum u Furijeovom spektru snage ukazuje na haotično ponašanje sistema. Vrednosti parametara su identične vrednostima na slici 7- 64d.\r\n8. ZAKLJUČAK\r\nIzvedeno istraživanje imalo je za cilj da utvrdi dinamička svojstva potresa izazvanih promenom naponskog stanja pri izradi horizontalnih podzemnih prostorija. Sastojalo se iz dva dela: u prvom delu izvršena je analiza vremenskih serija registrovanih potresa na površini terena, registrovanih u različitim rudnicima u Poljskoj (rudnici Rudna, Gornja Silezija i Legnica-Glogov). U drugom delu istraživanja, vršena je analiza dinamike fenomenološkog modela jednog i dva bloka sa oprugom.\r\nRezultati ispitivanja dinamike oscilovanja tla tokom rudarski generisanih potresa ukazali su na nelinearnu prirodu registrovanih serija, i to u pravcu sever-jug i u vertikalnom pravcu, dok oscilovanja tla u pravcu istok-zapad pripadaju klasi stohastičkih procesa sa Gausovom raodelom stohastičkog dela, koja može biti izmenjena nekom nepoznatom nelinearnom funkcijom. Međutim, i pored toga što oscilovanja tla u pravcu sever-jug i istok-zapad pripadaju grupi nelinearnih procesa, analiza nelinearnih vremenskih serija je pokazala da i u ovim pravcima oscilovanja tla pripadaju stohastičkim procesima. Ovi rezultati potvrđeni su malom vrednošću determinističkog faktora K (<1) i niskom greškom unakrsnog predviđanja kod testa stacionarnosti. Takođe, metoda prividno najbliže vrednosti nije dala rezultate u pogledu određivanja optimalne dimenzije razvijanja, što opet ukazuje na stohastičku prirodu registrovanih seizmograma.\r\nU drugom delu istraživanja, vršena je analiza modela bloka sa oprugom za različite sisteme, pri promeni vrednosti postojećih parametara, kao i uvođenjem novih parametara u model (jačina trenja c i vremensko kašnjenje T). U prvom slučaju, kod analize Karlson-Langer sistema, razmatrana je dinamika modela sa jednim i dva bloka. Analiza ovog modela izvedena je numerički, posmatranjem ponašanja rešenja sistema, uvodeći dva nova parametra: vremensko kašnjenje T i jačinu trenja c. Originalni model za jedan blok i fiksiranu vrednost bezdimenzionalne referentne brzine vc pokazuje samo oscilatorno ponašanje. Međutim, uvodeći prethodno pomenuta dva parametra, uočava se prelaz od ravnotežnog stanja ka periodičnom ponašanju. Takođe, uočeno je i da povećanjem vrednosti T, dolazi do prelaza modela iz periodičnog ponašanja u ravnotežno stanje (,,oscilatorna smrt''). Ova pojava je jako značajna, jer ukazuje na činjenicu da uvođenjem ovog parametra (vremensko kašnjenje T) ne dolazi samo do pojave kompleksnijih oblika dinamičkog ponašanja, već i do smirivanja bloka, odnosno njegovog kretanja ka ravnotežnom stanju. Pojava oscilatorne smrti takođe je osmatrana i za fiksiranu vrednost T, samo povećavajući parametar jačine trenja c.\r\nSlično dinamičko ponašanje se može uočiti i kod modela sa dva bloka. Originalni sistem, bez uvedenih parametara jačine trenja i vremenskog kašnjenja, takođe pokazuje samo oscilatorno ponašanje. Međutim, uvodeći prethodno pomenuta dva parametra, dolazi do prelaza iz periodičnog u kvaziperiodično ponašanje (i obrnuto) za različite vrednosti kontrolnih parametara. Ove oblasti različitog dinamičkog ponašanja su naročito dobro izražene u slučaju kada je prvi blok u ravnotežnom stanju, a parametri drugog bloka se menjaju. U ovom slučaju, ponašanje nalik determinističkom haosu osmatra se kao tranzijentna, netipična pojava.\r\nPrikazano istraživanje, iako sprovedeno samo numeričkim putem, ukazuje na nove pojave u proučavanju potresa izazvanih napredovanjem rudarskih radova. Prvo, uvodeći vremensko kašnjenje u izrazu za trenje, dinamičko ponašanje oba modela postaje kompleksnije. Izvedena analiza pokazuje da postoje dve jasno izražene i odvojene oblasti sa različitim dinamičkim ponašanjem - prva oblast za T<5, i druga oblast za T>8. S druge strane, uvođenjem parametra jačine trenja c, dolazi do pojave bifurkacije za relativno malu vrednost c (c <1). Ukoliko usvojimo da su parametar c i koeficijent trenja p. analogne veličine, rezultati izvedene analize bi mogli biti validni za plitke rasedne zone, na osnovu prethodno dobijenih rezultata u radu [Mizoguchi, i dr., 2007] gde je ustanovljena relativno niska vrednost koeficijenta trenja (0,50-0,58) duž Atotsugawa raseda u centralnom Japanu, na dubinama u rasponu 1-9km. Takođe, niske vrednosti parametra jačine trenja c ukazuju na relativno malu debljinu rasedne zone, što odgovara pretpostavci u radu [Marone, i dr., 1990], koji su pokazali da koeficijent trenja raste sa debljinom zdrobljenog rasednog materijala, uzrokujući promenu koeficijenta trenja p. i za više od dva puta. U pogledu mineralnog sastava duž rasedne zone, rezultati izvedene analize odgovaraju zdrobljenom kataklastičnom rasednom materijalu, izgrađenom od minerala sa relativno niskim koeficijentom trenja (grafit, talk, montmorilonit, hrizotil) i sa redukovanom tvrdinom u zasićenom stanju (antigorit, lizardit, kaolinit, muskovit, hlorit), prema [Morrow, i dr., 2000; Behnsen, Faulkner, 2012).\r\nU drugom, Bekerovom sistemu, vršena je analiza modela posmatranjem ponašanja rešenja, i to pretpostavljajući da se parametri menjaju periodično. Do pojave determinističkog haosa dolazi pri spregnutoj periodičnoj perturbaciji dva parametra - krutosti opruge, i promene čvrstoće tokom kretanja bloka. Neophodno je naglasiti da u ovom slučaju do pojave determinističkog haosa dolazi za mnogo veće vrednosti krutosti opruge, nego u radu [Becker, 2000].\r\nU trećem slučaju, vršena je analiza Madariaga sistema, i to na dva načina. Prvo, pretpostavljeno je da određeni parametri sistema (e i menjaju svoje vrednosti periodično (sinusoidalno), što se može pripisati dejstvu nekog spoljašnjeg faktora (vibracije mašina, efekti miniranja, i sl.), pri napredovanju rudarskih radova. Izvedena analiza je pokazala da promenom ugaone frekvencije oscilacija samo jednog od parametara, e ili £ dolazi do pojave kompleksne dinamike, sa prelazom od periodičnog i kvaziperiodičnog ponašanja ka determinističkom haosu. Međutim, pojava ovako kompleksne dinamike u slučaju perturbacije samo jednog parametra uslovljena je relativno velikim amplitudama oscilacija, SE i % koje se realno ne mogu očekivati u slučaju napredovanja rudarskih radova i iskopa podzemnih prostorija. Imajući to u vidu, sprovedena je dodatna analiza, u okviru koje je vršeno istovremeno variranje amplituda oscilacija oba parametra e i £ Rezultati ove dodatne analize ukazuju na to da je pojava haotičnog ponašanja, za spregnuto dejstvo oscilacija dva parametra, moguća za manje vrednosti amplituda oscilacija. Drugim rečima, moguće je da vrlo male oscilacije pojedinih parametara, čije se frekvencije poklapaju sa frekvencijom aseizmičkog kretanja duž raseda, izazovu reaktiviranje kretanja duž postojećih raseda. Realno je očekivati da se ovakve oscilacije mogu javiti pri različitim vrstama rudarskih aktivnosti (različiti izvori veštačkih vibracija - miniranje, rad mašina, i sl.).\r\nDrugi deo istraživanja obuhvatio je bifurkacionu analizu sistema sa uvedenim vremenskim kašnjenjem u izrazu za trenje, kojim se kvalitativno modeluje uticaj napredovanja rudarskih radova. Naime, utvrđeno je da povećanjem vrednosti parametra T (vremenskog kašnjenja) može doći i do pojave kompleksnijeg ponašanja, ali i do ,,smirivanja'' kretanja bloka, nakon čega se blok opet nalazi u ravnotežnom stanju (,,oscilatorna smrt''). U slučaju kada dolazi do pojave kompleksnijeg ponašanja, utvrđeno je da kretanje bloka prolazi kroz niz bifurkacija, od ravnotežnog stanja preko periodičnog i kvaziperiodičnog kretanja do pojave determinističkog haosa (Ruel-\r\nTakens-Njuhaus put u haos). Koristeći Rušeovu teoremu, određen je uslov (odnos između parametara sistema) za koji je sistem stabilan za sve vrednosti uvedenog vremenskog kašnjenja T.\r\nRezultati izvedenog istraživanja potvrđeni su činjenicom da je slično dinamičko ponašanje utvrđeno u modelima sa većim brojem blokova. Tako su, na primer, Erikson, i dr. [Erickson i dr., 2011] utvrdili postojanje determinističkog haosa za slične vrednosti kontrolnih parametara (s=0.5, = 0.5 and Y=0,5), ali u modelu sa više od 20 blokova. Nešto ranije istraživanje od strane Vieire [1999] ukazalo je na pojavu haosa u simetričnom modelu sa dva i tri bloka, ali za mnogo manje vrednosti kontrolnog parametra (1/vc ~ 0,122, gde je vc karakteristična brzina bloka).\r\nShodno tome, ovim istraživanjem je takođe potvrđeno da pojava determinističkog haosa u BK modelu sa Diterih-Ruina zakonom trenja nije zavisna od veličine modela, što je nije u saglasnosti sa rezultatima prethodnih istraživanja. Tako su, na primer, Schmittbuhl i dr. [Schmittbuhl, i dr., 1993] takođe izučavali model bloka sa oprugom i utvrdili da pojava haosa zavisi od veličine sistema. Slično, Erikson, i dr. [Erickson i dr. 2011] su utvrdili da kritična veličina parametra e potrebna za pojavu haosa opada kao funkcija broja blokova. Objašnjenje ovih rezultata može ležati u činjenici da samo zanemarivanje efekta ,,kašnjenja'', koje predstavlja inherentno svojstvo ne samo rudarski generisanih potresa, već i tektonskih zemljotresa, uzrokuje nemogućnost pojave kompleksnijeg ponašanja za sistem sa jednim blokom.\r\nNeophodno je naglasiti da je modelovanje tako složenih sistema samo jednim blokom sa oprugom u velikoj meri kvalitativno, i moguće je, u tom slučaju, možda da se pojave i efekti koji inače nisu prisutni in situ. Međutim, istraživanja modela sa većim brojem blokova, samo potvrđuju ispravnost našeg postupka. Na taj način, ostvaren je primarni cilj modelovanja - stvoriti fizički što jednostavniji model, a zadržati ponašanje sistema u što većoj mogućoj meri.\r\nU pogledu značaja izvedene analize za rudarski generisane potrese, doprinos ovih istraživanja ogleda se u sledećem:\r\n- analizom vremenskih serija registrovanih potresa na tri različite lokacije rudnika, potvrđena je stohastička priroda oscilovanja tla. Na taj način, opravdana je upotreba stohastičkih empirijskih obrazaca radi procene brzine i ubrzanja oscilovanja tla tokom potresa;\r\n- analiza dinamike sistema za model bloka sa oprugom, sa dva različita zakona trenja, i sa uvedenim parametrom jačine trenja c i vremenskim kašnjenjem T, potvrdila je kompleksnost dinamike modela, kroz niz bifurkacija, od ravnotežnog stanja, preko periodičnog i kvaziperiodičnog ponašanja do pojave determinističkog haosa. Neophodno je naglasiti da ova razlika u ponašanju modela mehanizma nastanka potresa (deterministički haotično) i registrovanog ubrzanja oscilovanja tla na površini terena (stohastičko) proizilazi iz sledećeg:\r\nO samo deo energije generisane kretanjem duž raseda biva emitovan u vidu seizmičkih talasa. Značajan deo energije ,,troši'' se na frikciono zagrevanje i topljenje, stvaranje silika-gela duž raseda, i dr. Na taj način, registrovana oscilovanja tla na površini ne oslikavaju u punoj meri kretanje duž raseda; O na svom putu do površine terena seizmički talasi gube energiju, i u zavisnosti od lokalne geološke građe terena, tektonike, debljine površinski raspadnute zone, zasićenosti tla vodom, topografije terena, kao i uticaja antropogenog faktora, seizmički talasi bivaju, u većoj ili manjoj meri, prigušeni. Kao takvi, mogu pokazivati stohastičko ponašanje.\r\n- u ovoj doktorskoj disertaciji definisana je metodologija izučavanja rudarski izazvanih potresa sa stanovišta nelinearne dinamike, koja obuhvata nekoliko faza: analizu nelinearnih vremenskih serija registrovanih potresa, numeričku analizu ponašanja rešenja sistema modela sa jednim i dva bloka u uslovima periodičnog oscilovanja parametara sistema, zatim modela sa uvedenim vremenskim kašnjenjem i parametrom jačine trenja, kao i bifurkacionu analizu modela mehanizma nastanka potresa sa uvedenim vremenskim kašnjenjem;\r\n- određene vrednosti parametara pri kojima dolazi do bifurkacija, odnosno do kvalitativne promene dinamike posmatranog sistema mogu biti važne, jer ukazuju na kritične vrednosti realno osmatranih veličina, koje su u ovom slučaju različite funkcije napona, a do čije promene neminovno dolazi napredovanjem rudarskih radova;\r\n- neophodno je napomenuti da metode i tehnike nelinearne dinamike i teorije haosa do sada nisu primenjivane u okviru ove naučne oblasti, što takođe predstavlja još jedan naučni doprinos u oblasti rudarskog inženjerstva. Rezultati izvedenih istraživanja daju podlogu daljim istraživanjima, u okviru kojih bi se u model bloka sa oprugom eksplicitno uveo parametar koji oslikava promenu naponskog stanja oko raseda. Usvajanjem funkcije njegove promene (linearne ili nelinearne) i to za one vrednosti koje su realno moguće, može se utvrditi veličina promene napona neophodna za generisanje rudarskih potresa.\r\nPRILOG 1. Rudarske oblasti u svetu u kojima su registrovani potresi izazvani različitom rudarskom aktivnošću.\r\nPRILOG 2. Algoritam metodologije nelinearnog dinamičkog modelovanja potresa izazvanih promenom naponskog stanja pri izradi horizontalnih podzemnih prostorija.\r\nPRILOG 3. Rezultati analize vremenske serije registrovanog potresa u basenu uglja u Gornjoj Sileziji u Poljskoj.\r\nSlika 1. Vremenske serije tipičnog ubrzanja oscilovanja tla za vreme potresa magnitude M~1,5 (£=1x107J) izazvanog rudarskim aktivnostima u basenu uglja u Gornjoj Sileziji, u tri pravca: (a) - sever-jug (748); (b) - istok zapad (865) i (c) - vertikalno (1126) [Dulinska, Fabijanska, 2011]\r\nSlika 2. Furijeov spektar snage za ubrzanje tla za vreme potresa magnitude M~1,5 (£=1x107J) izazvanog rudarskim aktivnostima u basenu uglja u Gornjoj Sileziji u (a) horizontalnom (sever-jug) i (b) vertikalnom pravcu. Kvalitativno sličan spektar snage se dobija i za ubrzanje tla u horizontalnom pravcu istok-zapad. Kontinualni spektar snage u oba slučaja ukazuje ili na moguće deterministički haotično ili na nasumično ponašanje posmatranog sistema.\r\nANALIZA SUROGAT PODATAKA\r\nSlika 3. Testiranje I nulte hipoteze: (a) pravac sever-jug; (b) pravac istok-zapad; (c) vertikalno. Očigledno je da u sva tri slučaja nulta hipoteza može biti odbačena, s obzirom na to da je e0