Katarzyna SZAFULERA, Marek JENDRYŚ Politechnika Śląska, Gliwice
Katedra Geomechaniki, Budownictwa Podziemnego i Zarządzania Ochroną Powierzchni
NUMERYCZNA ANALIZA STANU ODKSZTAŁCENIA GÓROTWORU W SĄSIEDZTWIE PŁYTKICH WYROBISK PORUDNYCH
Streszczenie. W artykule została przedstawiona numeryczna analiza stanu odkształcenia górotworu w sąsiedztwie płytkich wyrobisk porudnych. Deformacje górotworu w otoczeniu tych wyrobisk zostały określone przy użyciu programu FLAC3D. W wyniku modelowania numerycznego zostały określone wartości poszczególnych składowych tensora odkształcenia, a w konsekwencji wartości i kierunki odkształceń głównych.
NUMERICAL ANALYSIS OF ROCK MASS STRAIN STATE IN THE VICINITY OF SHALLOW ORE EXCAVATIONS
Summary. Paper presents numerical analysis of rock mass strain state in the vicinity of shallow ore excavations. Deformation of rock mass has been calculated by use of FLAC3D. As a result of numerical modeling the value of strain tensor components has been identified and consequently the values and directions of principal strain has been determined.
1. Wprowadzenie
Pustki pogórnicze sprzyjające powstawaniu deformacji nieciągłych zagrażają bezpieczeństwu użytkowania powierzchni oraz stanowią zagrożenie dla środowiska przyrodniczego.
Szczególnie zagrożone występowaniem deformacji nieciągłych w wyniku reaktywacji pustek po płytkiej eksploatacji rudnej są rejony Bytomia, Tarnowskich Gór, a także Piekar Śląskich, należące, zgodnie z pracą [8], do obszaru dużego zagrożenia zapadliskowego – kategoria D.
K. Szafulera, M. Jendryś 188
Opierając się na udokumentowanych deformacjach nieciągłych, wielu badaczy − m.in.
M. Chudek i in. (1988) oraz W. Janusz i A. Jarosz (1976) − wyraziło pogląd, że warunkiem koniecznym powstania powierzchniowej deformacji nieciągłej jest samopodsadzenie pustki istniejącej w górotworze. Może ono powstać w sposób gwałtowny, a czas jego wystąpienia jest trudny do określenia i może być przesunięty w odniesieniu do czasu wytworzenia się pustki nawet o kilkadziesiąt lat.
Według M. Boreckiego (1980) i M. Chudka (1981) jedną z przyczyn utarty stateczności wyrobisk są zmiany stanów naprężeniowo-deformacyjnych, wywołane eksploatacją górniczą, przy czym przekroczenie granicznych wartości odkształceń powoduje zniszczenie struktury skał i wytworzenie w otoczeniu wyrobiska strefy spękań. Skały spękane pod wpływem ciężaru własnego mogą ulec przemieszczeniu w kierunku wyrobisk i spowodować ich samopodsadzenie. W świetle powyższego można przyjąć, że występowanie deformacji nieciągłych na powierzchni jest uzależnione od stateczności płytko zalegających wyrobisk, dla których zasadniczym zagadnieniem jest określenie stanów naprężeniowo-odkształceniowych górotworu w otoczeniu analizowanych wyrobisk.
W artykule przedstawiono, uzyskane na podstawie modelowania numerycznego, wyniki obrazujące stan odkształcenia górotworu w otoczeniu wyrobisk zlokalizowanych na małych głębokościach, nieprzekraczających 100 m.
2. Modelowanie numeryczne górotworu
Do zbadania odkształceń wywołanych płytką eksploatacją rud cynku i ołowiu został wybrany obszar górotworu, w rejonie którego w okresie kilkunastu lat zinwentaryzowano na powierzchni 7 deformacji nieciągłych (zapadlisk).
2.1. Przestrzenny model górotworu
Dla przeprowadzenia symulacji numerycznych został zbudowany przestrzenny model, tworzący prostopadłościenny blok górotworu o wysokości równej 100 m, podstawie w kształcie kwadratu o boku równym 500 m (rys. 1). Objętość zamodelowanej bryły górotworu jest równa 25 000 000 m3 .
500 500
100
LEGENDA:
- glina piaszczysta - glina
- piasek
dolomit kremowo -
- ił witriolowy żółty silnie zw.
-
dolomit czerwony zwietrzały - dolomit brunatnoczerwony
- rudy utlenione rudy siarczkowe - wapień krystaliczny
Rys. 1. Widok ogólny modelowanej bryły górotworu Fig. 1. General view of analyzed rock mass
Zarówno model górotworu, jak i przestrzenny model warstwy złoża zostały zbudowane z predefiniowanych, prostopadłościennych stref o nazwie brick, dostępnych w programie FLAC3D.
Wymiary elementów w poszczególnych rejonach zostały dobrane tak, aby ich największe zagęszczenie występowało w tych obszarach górotworu, gdzie znajdują się wyrobiska, nad którymi powstały trzy rozpatrywane zapadliska.
Dane geologiczne dotyczące litologii i miąższości warstw skalnych w modelowanej bryle górotworu zostały ustalone na podstawie profili litologicznych szybów, zlokalizowanych najbliżej analizowanego obszaru (szyb „P”, szyb „K”) oraz profilu sporządzonego na podstawie dokumentacji geologicznej byłych Zakładów Górniczo-Hutniczych „OB”.
2.2. Przestrzenny model warstwy złoża
Na podstawie archiwalnej mapy wyrobisk rudnych został opracowany przestrzenny model warstwy złoża, tworzący prostopadłościenny blok o wysokości 4,5 m, podstawie w kształcie kwadratu, którego bok równy jest 500 m, i objętości wynoszącej 1 125 000 m3. Eksploatację rudną zasymulowano w środkowej części warstwy złoża, w prostopadłościennym bloku o wymiarach 285 x 255 x 4,5 m.
K. Szafulera, M. Jendryś 190
Rys. 2. Widok archiwalnej mapy wyrobisk górniczych w warstwie rudy z podziałem na parcele Fig. 2. The view of archival map of mining excavations in the layer of ore divided into parcels
W przedmiotowym rejonie rudy cynku i ołowiu wybierano systemem filarowo- komorowym. Polegał on na drążeniu z chodnika eksploatacyjnego wąskich zabierek o maksymalnej szerokości 5 m i długości do 40 m, przy stałej grubości wybieranej warstwy równej 4,5 m. Była to eksploatacja jedno- i dwuskrzydłowa, najczęściej prowadzona diagonalnie.
W celu dokładnego odwzorowania skomplikowanej sytuacji górniczej warstwa złoża została podzielona na wiele parcel (rys. 2). Znając miąższość rudy oraz współrzędne parcel wyznaczono prostopadłościenne regiony zbudowane z predefiniowanych w programie FLAC3D stref o nazwie brick. Na rys. 3 przedstawiono kształty i wielkości regionów tworzących warstwę złoża, a uzyskaną po zakończeniu eksploatacji sytuację górniczą przedstawia rys. 4.
X Y
Z
Block Group
r0 r1 r2 r3 r4 r5 r6 r7 r8 r9 r10 r11 r12 r13 r14 r15 r16 r16' r17 r18
Rys. 3. Widok warstwy złoża z podziałem na regiony Fig. 3. The view of ore deposit with an indication of regions
Block Group
niepodsadzone podsadzka_piaskowa zroby_zawaŁowe ruda
Rys. 4. Widok warstwy złoża po zakończeniu eksploatacji Fig. 4. The view of ore deposit after exploitation
2.3. Przyjęte modele konstytutywne
W modelowanej bryle górotworu przyjęto kilka modeli fizykalnych. Występującym przy powierzchni gruntom został przyporządkowany model ciała sprężysto-plastycznego ze zmodyfikowanym warunkiem wytrzymałościowym Coulomba–Mohra [6]. Dla triasowego górotworu zwięzłego został przyjęty model sprężysto–plastyczny z osłabieniem materiału skalnego w stadium pokrytycznym (strain-softening) [6]. Rumowisko zawałowe zostało opisane modelem sprężystym, natomiast pustkom poeksploatacyjnym wypełnionym podsadzką piaskową został przyporządkowany model sprężysto-plastyczny ze wzmocnieniem w stadium pokrytycznym (strain-ardening) [6].
2.4. Stałe materiałowe modelowanego górotworu
Badania skał w rejonie płytkiej eksploatacji rudnej były prowadzone przez Kuhla w 1954 roku i obejmowały wyrobiska usytuowane w północno–środkowej części Bytomia oraz w dzielnicy Szarlej miasta Piekary Śląskie. W 1965 roku Stacja Badawczo- Doświadczalna Górnictwa Rud przy ZGH „Orzeł Biały” przeprowadziła podobne badania w wyrobiskach zlokalizowanych w środkowej i wschodniej części Bytomia. Otrzymane z tych dwóch źródeł wyniki są porównywalne.
Dysponując wynikami badań własności fizykomechanicznych oraz szczegółowym opisem budowy geologicznej, dla właściwej oceny jakości analizowanego górotworu, czyli dopasowania wartości stałych materiałowych dla skał w modelowanej numerycznie bryle
K. Szafulera, M. Jendryś 192
górotworu, zastosowano klasyfikację E. Hoeka – Geological Strength Index [4]. Wartości stałych materiałowych dla sześciu różnych rodzajów skał, ustalonych na podstawie klasyfikacji GSI, zostały przedstawione w tablicy 1. Stałe materiałowe dla gruntów zostały ustalone na podstawie danych opublikowanych w pracach [9] oraz [10].
Tablica 1 Wprowadzone do modelu numerycznego wartości stałych materiałowych gruntów i skał, występujących w analizowanym obszarze
3. Komputerowa symulacja eksploatacji filarowo-komorowej
Harmonogram symulacji eksploatacji rud cynku i ołowiu w obrębie modelowanej bryły górotworu został ustalony przy wykorzystaniu archiwalnej mapy wyrobisk górniczych.
Biorąc pod uwagę rodzaj symulowanych wyrobisk eksploatacyjnych, założono, że w jednym stadium symulacji ruda zostaje wybrana, uruchamiane są obliczenia, w czasie których model zostaje doprowadzony do stanu równowagi i dopiero wówczas, w kolejnym stadium symulacji następuje wypełnienie pustki poeksploatacyjnej.
Efekt konsolidacji gruzowiska skalnego w strefie zawału oraz przyrost gęstości i zmiany stałych sprężystości podsadzki piaskowej w czasie były symulowane w kilku stadiach: dla zawału – 4 stadia, dla podsadzki – 3 stadia.
W podstawowej wersji programu FLAC3D istnieje możliwość odczytywania uzyskanych wyników obliczeń m.in. w postaci: składowych stanu naprężenia, przemieszczeń liniowych, wektorów przemieszczeń, przyrostów i prędkości odkształceń postaciowych oraz odkształceń objętościowych. Brak natomiast informacji na temat składowych tensora odkształcenia i odkształceń głównych.
Zgodnie ze wskazówkami twórców programu [5], w wewnętrznym języku programowania został utworzony algorytm umożliwiający prezentację i zapis wszystkich składowych tensora odkształcenia
( )
Tε dla dowolnie wybranego obszaru modelu. Dla zdefiniowanego przez użytkownika punktu, określonego współrzędnymi (x,y,z), algorytm ten w programie pozwala odczytać składowe tensora odkształcenia, wyznaczone dla strefy różnic skończonych (ipz), której środek znajduje się najbliżej zdefiniowanego punktu.Zastosowanie tak opracowanego algorytmu pozwoliło poznać wartości odkształceń liniowych i postaciowych w modelowanej bryle górotworu. Na podstawie tych wartości, korzystając z równania kubicznego i niezmienników odkształcenia [12], [13], zostały określone największe (ε1), pośrednie (ε2) i najmniejsze (ε3) odkształcenia główne.
Tensor odkształcenia zostały odczytany z programu FLAC3D dla 7830 punktów. Punkty były węzłami przestrzennej regularnej siatki o wymiarach 285 m × 255 m. Skok po współrzędnych x, y był równy 10 m, a po współrzędnej z: 10 m, 5 m, 3 m i 2 m w zależności od głębokości (z = 100, 90, 80, 70, 60, 57, 52, 47, 45, 43 m). Widok przestrzennej siatki punktów przedstawiony został na rys. 7.
Otrzymane wyniki obliczeń w postaci map największego(ε1) i najmniejszego (ε3) odkształcenia głównego górotworu po zakończonej eksploatacji rud cynku i ołowiu zostały przedstawione na rys. 5 − 6 w płaszczyźnie pionowej, przebiegającej od warstwy wybranej rudy do powierzchni, oraz w płaszczyźnie stropu eksploatowanego złoża rudy (rys. 7, 8).
Wartości odkształceń są wyrażone w promilach, znakiem ujemnym zostały opatrzone skrócenia, znakiem dodatnim − wydłużenia.
K. Szafulera, M. Jendryś 194
Rys. 5. Rozkład największego odkształcenia głównego (ε1) w płaszczyźnie XZ, y = 48 m, modelowanej bryły górotworu, od stropu wyrobisk porudnych (z = 43 m) do powierzchni (z = 100 m) po zakończeniu eksploatacji
Fig. 5. Distribution of maximum principal strain (ε1) in the plain XZ, y = 48, in the vertical section from excavation roof to ground surface
Rys. 6. Rozkład najmniejszego odkształcenia głównego (ε3) w płaszczyźnie XZ, y = 48 m, modelowanej bryły górotworu, od stropu wyrobisk porudnych (z = 43 m) do powierzchni (z = 100 m), po zakończeniu eksploatacji
Fig. 6. Distribution of minimum principal strain (ε3) in the plain XZ, y = 48, in the vertical section from excavation roof to ground surface
Rys. 7. Rozkład największego odkształcenia głównego (ε1) w płaszczyźnie stropu złóż rud cynku i ołowiu w płaszczyźnie XY, z = 43 m, modelowanej bryły górotworu po zakończeniu eksploatacji
Fig. 7. Distribution of maximum principal strain (ε1) in the plain of roof of ore deposit
Rys. 8. Rozkład najmniejszego odkształcenia głównego (ε3)w płaszczyźnie stropu złóż rud cynku i ołowiu w płaszczyźnie XY, z = 43 m, modelowanej bryły górotworu po zakończeni eksploatacji
Fig. 8. Distribution of minimum principal strain (ε3) in the plain of roof of ore deposit
K. Szafulera, M. Jendryś 196
4.1. Analiza uzyskanych numerycznie wyników obliczeń
W wyniku prowadzonej w latach 1945 – 1984 w przedmiotowym rejonie eksploatacji rud cynku i ołowiu została wykonana, a następnie zlikwidowana lub pozostawiona bez podsadzenia bardzo duża liczba wyrobisk o zróżnicowanych wymiarach. Ta niezwykle złożona sytuacja górniczo-geologiczna była przyczyną deformacji skał nadległych, a szczególnie warstw zalegających bezpośrednio nad prowadzoną eksploatacją.
W częściach zachodnio-centralnej i południowej analizowanej bryły górotworu wystąpiły maksymalne odkształcenia. W miejscach tych największe skrócenia osiągnęły wartości 0,66 mm/m. Wydłużenia wystąpiły w dwóch kierunkach, a ich największa wartość wynosiła 0,26 mm/m. Deformacje te były skutkiem prowadzonej w latach 1945 – 1956 eksploatacji rud z zawałem stropu. Strefy dużych deformacji wystąpiły nad zrobami zawałowymi (rys. 5, 6) i objęły znaczną część górotworu zarówno w płaszczyźnie poziomej, jak i pionowej. Strefy największych skróceń sięgnęły do ok. 14 m nad stropami wyrobisk, a powyżej tej wysokości malały, osiągając na powierzchni wartość ok. 0,21 mm/m (rys. 6). Największe wydłużenia objęły swym zasięgiem partię górotworu do ok. 20 m nad stropami wyrobisk (0,26 mm/m) i kolejno malały w kierunku powierzchni do wartości 0,13 mm/m (rys. 5).
Równie duże odkształcenia wystąpiły tam, gdzie w latach 1972 – 1984 była prowadzona zgrupowana eksploatacja z ochroną stropu, jednak zasięgi tych deformacji nie były rozległe, miały charakter lokalny (rys. 7, 8).
W stropach wyrobisk, gdzie ruda była wybierana pojedynczymi zabierkami, podobnie w stropach niepodsadzonych chodników, odkształcenia w dwóch kierunkach były skróceniami, a ich największa wartość wynosiła 0,35 mm/m (rys. 8). W stropach analizowanych wyrobisk największe wydłużenia osięgnęły wartość ok. 0,19 mm/m (rys. 7).
5. Podsumowanie
Jedną z metodą pozyskiwania informacji o procesach zachodzących na powierzchni i w masywie skalnym na skutek prowadzonych robót górniczych jest budowa i obserwacja modeli górotworu.
W celu odzwierciedlenia skomplikowanej, przestrzennej sytuacji górniczej wyrobisk o prostokątnym przekroju poprzecznym i obliczenia odkształceń górotworu, wywołanych
systemem filarowo-komorowym z zawałem skał stropowych oraz z ochroną stropu, prowadzona w latach 1945 – 1984.
Rezultatem tych symulacji są wyniki obliczeń w postaci pełnego tensora odkształceń górotworu w sąsiedztwie płytkich wyrobisk porudnych oraz wyznaczone na podstawie podanych w pracy wzorów odkształcenia główne górotworu.
Prawidłowo przeprowadzona numeryczna analiza stanu odkształcenia górotworu może stanowić jeden z punktów wyjściowych w badaniu stateczności górotworu w sąsiedztwie płytkich (zwłaszcza dawnych) wyrobisk powęglowych i porudnych i występowania deformacji nieciągłych.
BIBLIOGRAFIA
1. Borecki M.: Mechanika budowli podziemnych. Skrypt Politechniki Śląskiej, nr 904, Gliwice 1980.
2. Chudek M., Janusz W., Zych J.: Studium dotyczące stanu rozpoznania tworzenia się i prognozowania deformacji nieciągłych pod wpływem podziemnej eksploatacji złóż. Zeszyty Naukowe Politechniki Śląskiej, s. Górnictwo, z. 141, Gliwice 1988.
3. Chudek M.: Mechanika górotworu. Skrypt Politechniki Śląskiej, nr 956/61, Gliwice 1981.
4. Hoek E., Brown E. T.: Practical estimates of rock mass strength. Int. J. Rock Mechanics and Min. Sci. and Geomech. Abst., Vol. 34, No. 8, 1988.
5. Itasca Consulting Group Inc., 3DEC, User’s Guide, Minneapolis b 2005.
6. Itasca Consulting Group Inc., FLAC 3D, User’s Guide, Minneapolis a 2005.
7. Janusz W., Jarosz A.: Nieciągłe deformacje powierzchni terenu wywołane płytką podziemną eksploatacją górniczą. Państwowa Akademia Nauk, Katowice 1976
8. Kotyrba A., Tyrała A., Michalak J.: Płytka eksploatacja złóż węgla kamiennego i rud cynkowo-ołowiowych oraz metody eliminacji jej wpływu na warunki geologiczno- inżynierskie województwa katowickiego. Prace Naukowe Głównego Instytutu Górnictwa, s. Konferencje, nr 24, 1998.
9. Kwaśniewski M.: Ocena stateczności masywu skalnego w sąsiedztwie uskoków poddanych wpływom eksploatacji górniczej (na przykładzie uskoku III i IV w KWK
„Pokój” i „Polska Wirek”, Gliwice 1996 (praca niepublikowana).
10. Smolnik G.: Właściwości procesu osiadania terenu górniczego w czasie. Praca doktorska (niepublikowana). Politechnika Śląska, Wydział Górnictwa i Geologii, Gliwice 2009.
11. Szafulera K.: Wpływ eksploatacji górniczej na stateczność płytkich wyrobisk porudnych i występowanie deformacji nieciągłych. Praca doktorska (niepublikowana). Politechnika Śląska, Wydział Górnictwa i Geologii, Gliwice 2011.
12. Walczak J.: Wytrzymałość materiałów oraz podstawy teorii sprężystości i plastyczności.
Tom I. Państwowe Wydawnictwo Naukowe, Łódź – Warszawa – Kraków 1963.
K. Szafulera, M. Jendryś 198
13. Walczak J.: Wytrzymałość materiałów oraz podstawy teorii sprężystości i plastyczności.
Tom II. Państwowe Wydawnictwo Naukowe, Łódź – Warszawa – Kraków 1959.
Abstract
Paper presents numerical analysis of rock mass strain state in the vicinity of shallow ore excavations. Deformation of rock mass has been calculated by use of FLAC3D. Numerical model used in simulation reflecting the state of rock mass in the vicinity of shallow excavation and has been build on the base archival mining maps and litological profiles of shaft localized in analyzed area.
As results of numerical modeling the value of strain tensor components has been identified and consequently the values and directions of principal strain have been determined.
The results of numerical calculation have been presented in the form of distributions of major and minor principal strain. These distributions of strain state have been presented for chosen section of numerical model.
Obtained results are input data for further effort analysis and prediction of rock mass behavior, especially for assessment of occurrence of surface deformation.
