Index of /rozprawy2/10701

Pełen tekst

(1)AKADEMIA GÓRNICZO-HUTNICZA IM. STANISŁAWA STASZICA WYDZIAŁ WIERTNICTWA NAFTY I GAZU KATEDRA INŻYNIERII GAZOWNICZEJ. ROZPRAWA DOKTORSKA. WPŁYW PARAMETRÓW TERMODYNAMICZNYCH PROCESU MAGAZYNOWANIA GAZU NA POJEMNOŚĆ KAWERNY W ZŁOŻU SOLI KAMIENNEJ. MGR INŻ. KAROLINA SERBIN. PROMOTOR: DR HAB. INŻ. JAROSŁAW ŚLIZOWSKI, PROF. AGH. KRAKÓW 2013.

(2) Spis treści SPIS WAŻNIEJSZYCH OZNACZEŃ ............................................................................................3 1. WPROWADZENIE ................................................................................................................5 1.1. Proces magazynowania gazu w kawernach solnych......................................................5 1.2. Parametry termodynamiczne mające wpływ na pojemność kawerny ............................6 1.3. Aktualny stan wiedzy...................................................................................................7 1.4. Cel i zakres pracy.........................................................................................................9 2. TEORETYCZNE PODSTAWY MODELI TERMODYNAMICZNYCH KAWERN MAGAZYNOWYCH GAZU .......................................................................................................12. 2.1. Procesy zachodzące w magazynowanym gazie...........................................................13 2.1.1. Równanie stanu magazynowanego gazu..............................................................13 2.1.2. Równanie ciągłości gazu .....................................................................................15 2.1.3. Równanie ruchu gazu ..........................................................................................15 2.1.4. Równanie energii gazu ........................................................................................17 2.1.5. Równanie oporu przepływu przez kolumnę rur....................................................19 2.2. Wymiana ciepła pomiędzy gazem a górotworem........................................................20 2.3. Zmiany temperatury w górotworze solnym – równanie przewodnictwa cieplnego......21 2.4. Algorytm programu KAGA .......................................................................................22 2.4.1. Opis segmentu zatłaczania gazu RUZA...............................................................22 2.4.2. Opis segmentu PRZERWA .................................................................................25 2.4.3. Opis segmentu poboru gazu z kawerny RUPO ....................................................27 3. TEORETYCZNE PODSTAWY MODELI GEOMECHANICZNYCH KAWERN MAGAZYNOWYCH GAZU .......................................................................................................29. 3.1. Równanie konstytutywne soli kamiennej....................................................................29 3.1.1. Odkształcenia sprężyste ......................................................................................31 3.1.2. Odkształcenia pełzania........................................................................................32 3.2. Temperatura i naprężenie pierwotne w górotworze ....................................................33 3.3. Algorytm programu GEOSOLK ................................................................................34 4. LABORATORYJNE BADANIA PRĘDKOŚCI PEŁZANIA SOLI ...................................................38. 1.

(3) 5. ZAKRES OBLICZEŃ NUMERYCZNYCH ................................................................................44 5.1. Warunki brzegowe i dyskretyzacja obszaru dla programu KAGA ..............................44 5.2. Warunki brzegowe i dyskretyzacja obszaru dla programu GEOSOLK .......................46 5.3. Rozpatrywane prawa pełzania....................................................................................46 5.4. Wychłodzenie górotworu podczas ługowania i scenariusz eksploatacji kawerny magazynowej....................................................................................................................49 5.5. Oznaczenie modeli numerycznych .............................................................................54 6. WPŁYW PARAMETRÓW TERMODYNAMICZNYCH PROCESU MAGAZYNOWANIA NA OBJĘTOŚĆ KAWERNY............................................................................................................56. 6.1. Analiza konwergencji i prędkości konwergencji kawerny ..........................................56 6.1.1. Wpływ wychłodzenia górotworu i prawa pełzania...............................................56 6.1.2. Wpływ czynnika temperaturowego Q/R ..............................................................64 6.1.3. Wpływ głębokości...............................................................................................71 6.2. Analiza przemieszczeń...............................................................................................74 7. PODSUMOWANIE ..............................................................................................................80 LITERATURA ........................................................................................................................83 SPIS RYSUNKÓW ...................................................................................................................86 SPIS TABEL ...........................................................................................................................89. 2.

(4) Spis ważniejszych oznaczeń A – stały współczynnik w prawie pełzania Nortona a, b – stałe, zależne od składu gazu Cg – ciepło właściwe górotworu D – średnica hydrauliczna rury w równaniu oporu przepływu przez kolumnę rur dl – elementy całkowania po długości rury dS – element całkowania po powierzchni w równaniu ciągłości gazu dV – element całkowania po objętości w równaniu ciągłości gazu dε – całkowity przyrost odkształceń. dε el – przyrost odkształceń sprężystych (elastycznych) dε v – przyrost odkształceń lepkich (pełzania) E – moduła Younga (odkształcalności podłużnej). e – energia właściwa gazu (całkowita) G=. E - moduł Kirchhoffa (odkształcalności postaciowej) 2(1 + ν ). g – wektor przyspieszenia ziemskiego |g|=g h – głębokość od powierzchni terenu K=. E - moduł ściśliwości 3(1 − 2ν ). k – szorstkość bezwzględna L – długość wycinka przepływu, tzn. rozpatrywanego odcinka rury w równaniu oporu przepływu przez kolumnę rur L – praca mechaniczna wykonana na gazie w równaniu energii gazu n – wektor jednostkowy, normalny do powierzchni Γ, w równaniu ciągłości gazu n – współczynnik potęgi w prawie pełzania Nortona p – ciśnienie Q – energia aktywacji R – uniwersalna stała gazowa Re – liczba Reynoldsa T – temperatura, Tg – temperatura górotworu Tk – temperatura wewnątrz kawerny 3.

(5) t – czas u – zredukowana energia wewnętrzna gazu v – prędkość przepływu gazu v – objętość molowa Z – powierzchnia boczna rozpatrywanego odcinka rury w równaniu oporu przepływu przez kolumnę rur yi – molowy udział i-tego składnika w mieszaninie gazowej, α – współczynnik wymiany cieplnej Γ – powierzchnia zewnętrzna elementu Ω w równaniu ciągłości gazu. δ ij – delta Kroneckera ε ef – odkształcenie efektywne. ε ijel – składowe tensora odkształceń ε m = 13 (ε 11 + ε 22 + ε 33 ) – odkształcenie średnie $I – tensor jednostkowy. σˆ – tensor naprężeń w gazie. σ ef – naprężenie efektywne σ ij – składowe tensora naprężeń σ m = 13 (σ 11 + σ 22 + σ 33 ) – naprężenie średnie λ – bezwymiarowy współczynnik oporu λg – współczynnik przewodnictwa cieplnego w górotworze $ – tensor naprężeń stycznych Λ. ν – współczynnik Poissona ρ – gęstość gazu ρg – gęstość górotworu Φ – strumień cieplny Ω – dowolny element kawerny lub otworu, wypełniony gazem w równaniu ciągłości gazu. 4.

(6) 1. Wprowadzenie. 1.1. Proces magazynowania gazu w kawernach solnych. Z punktu widzenia ilości magazynowanego gazu kawerny magazynowe w złożach soli kamiennej są trzecim, po sczerpanych złożach gazu i ropy naftowej oraz szczelnych strukturach wodonośnych tzw. aquiferach, sposobem podziemnego magazynowania tego surowca. Zgodnie z (Reinisch, 2000) ilość aktywnego gazu zmagazynowana na świecie w poszczególnych typach magazynów przedstawia się następująco: sczerpane złoża 192,4 mld m3 (76,47%), struktury zawodnione 47,6 mld m3 (18,93%), kawerny solne 11,5 mld m3 (4,57%), kopalnie i inne obiekty 0,1 mld m3 (0,03%). Z punktu widzenia zaspokojenia niedoborów szczytowych rola kawern jest znacznie większa a głównymi ich zaletami w porównaniu z innymi rodzajami magazynów są: • znacznie większa i stała w czasie wydajność możliwego poboru i zatłaczania gazu, • możliwość wielu cykli napełniania i opróżniania w ciągu roku, • mała ilość gazu buforowego nie przekraczająca zwykle 20% objętości roboczej, • mniejszy koszt inwestycji w przeliczeniu na jednostkę wydajności poboru gazu z magazynu. Podstawowymi wadami są natomiast wyższe koszty całego magazynu na jednostkę objętości magazynowanego gazu i dłuższy czas trwania budowy. Kawerny wykonywane są z powierzchni tzw. metodą otworową polegającą na rozpuszczaniu soli kamiennej w złożu wodą i wydobywaniu jej w postaci solanki. Warunki eksploatacji komory magazynowej gazu są całkowicie odmienne od normalnej eksploatacji otworowej i eksploatacji komór magazynowych ciekłych węglowodorów gdzie ciśnienie w komorach ulega niewielkim wahaniom pomiędzy ciężarami słupa solanki i magazynowanego medium. Komory magazynowe gazu są zbiornikami oddychającymi. Po zakończeniu ługowania do komory zapuszczane są dwie gazoszczelne kolumny rur, pierwsza to sięgająca stropu komory tzw. kolumna wydobywcza gazu, a druga to sięgająca spągu komory kolumna solankowa do opróżniania komory z solanki, które jest zarazem jej pierwszym napełnianiem, po którym kolumna solankowa jest wyciągana, albo cięta i zrzucana na dno kawerny. Później magazyn funkcjonuje pod zmiennym ciśnieniem, którego wartości skrajne określa się indywidualnie dla każdej komory.. 5.

(7) Zmienność ciśnienia w kawernie powoduje zmienność temperatury magazynowanego gazu, a tym samym na skutek wymiany ciepła zmienność temperatury otaczającego górotworu. Efektem tego jest zmienna prędkość pełzania górotworu, która zależy zarówno od ciśnienia gazu w komorze jak i temperatury górotworu. Sól kamienna jest bowiem ośrodkiem o bardzo silnych właściwościach reologicznych, którego odkształcenia zmieniają się w czasie nawet przy niezmienionym poziomie naprężeń i temperatury. Formułę opisującą te zmiany nazywamy prawem pełzania, które wyraża prędkość odkształceń efektywnych jako funkcję naprężeń efektywnych, czasu lub osiągniętych odkształceń efektywnych, temperatury, wilgotności i parametrów strukturalnych. Problem wyboru odpowiedniej formuły prawa pełzania nie został rozwiązany jednoznacznie i w literaturze można spotkać przykłady znacznie różniące się stopniem komplikacji.. 1.2. Parametry termodynamiczne mające wpływ na pojemność kawerny Parametry termodynamiczne odgrywają podstawową rolę w procesie podziemnego magazynowania gazu ziemnego w kawernach solnych. Należy brać je pod uwagę zarówno podczas projektowania magazynu kawernowego, jak i podczas jego eksploatacji, aby wybrać optymalne szybkości opróżniania i napełniania poszczególnych kawern. Do termodynamicznych zjawisk mających główny wpływ na zachowanie się kawerny magazynowej należą: • wychłodzenie górotworu w trakcie ługowania kawerny spowodowane niższą temperaturą medium ługującego (woda lub solanka nienasycona) w porównaniu do temperatury pierwotnej górotworu, • termo - hydrodynamiczne przemiany zachodzące w gazie ziemnym zatłaczanym, przechowywanym lub pobieranym z kawerny, których efektem jest wzrost temperatury gazu w trakcie zatłaczania i jej spadek przy opróżnianiu komory, • wymiana ciepła pomiędzy gazem ziemnym a górotworem solnym otaczającym kawernę i otwór prowadzący do kawerny w trakcie eksploatacji kawerny, • zmiany temperatury w górotworze solnym związane z rozchodzeniem się ciepła na skutek przewodnictwa cieplnego górotworu.. 6.

(8) Bezpośredni pomiar tych procesów w kawernie i górotworze praktycznie nie jest możliwy ze względu na brak odpowiednich metod pomiaru. O ich przebiegu możemy wnioskować jedynie na podstawie danych pośrednich, takich jak dane rejestrowane podczas operacji magazynowych na głowicy otworu. Mierzy się tam trzy podstawowe parametry gazu: ciśnienie, temperaturę i wydajność przepływu. Wprawdzie teoretycznie rzecz biorąc można wewnątrz kawerny zainstalować urządzenia mierzące w sposób ciągły temperaturę i ciśnienie, ale wyprowadzenie połączeń rejestratorów na powierzchnię poprzez wgłębny zawór bezpieczeństwa bez zakłócania jego pracy jest skrajnie utrudnione i kosztowne. Nie ma natomiast możliwości mierzenia rozkładu temperatury w górotworze wokół kawerny, ani też zmierzenia in situ rzeczywistych współczynników wymiany ciepła przez jej powierzchnię. Jedynym sposobem, aby otrzymać informacje o stanie cieplnym kawerny i otaczającego ją solnego górotworu jest posłużenie się modelem komputerowym, w którym we właściwy sposób będą symulowane wszystkie zjawiska termo- i hydrodynamiczne związane z magazynowaniem gazu. Zasadniczym celem modelowania jest: uzyskanie rozkładu temperatur w górotworze, które wpływają na prędkość jego pełzania a w konsekwencji konwergencje kawerny.. 1.3. Aktualny stan wiedzy Termodynamiczne zjawiska zachodzące w gazie magazynowanym w kawernie oraz termomechaniczne w otaczającym ją górotworze są w zasadzie dobrze znane od kilkunastu lat, niemniej jako zjawiska odrębne, każde z osobna, nie były natomiast rozpatrywane jako wzajemnie wpływające na siebie procesy. Na przykład dobrze znane jest równanie przewodnictwa cieplnego, jak też fakt zależności współczynników przewodnictwa od stanu naprężeń, ale nikt nie opracował odpowiednich formuł ilościowych i nie analizował propagacji ciepła w górotworze solnym przy zmieniających się naprężeniach. Wiadomo również że temperatura ma bardzo istotny wpływ na prędkość pełzania stacjonarnego soli kamiennej. Również dobrze znane są równania stanu gazu ziemnego i jego zachowanie się w dużych zbiornikach, czy przy przepływie przez rury, długi czas jednak mało było wiadomo o wymianie ciepła między magazynowanym gazem i górotworem oraz jak wpływa ona na stan gazu w kawernie. Dopiero prace Thaulego (Thaule, 1994), (Thaule, 1997) i stworzenie. 7.

(9) programu Kameleon II umożliwiającego szczegółowe modelowania na gęstych siatkach aproksymacyjnych obejmujących warstwę przyścienną z możliwością uwzględnienia kilku modeli turbulencji, wyjaśniły które wymiary powinny być stosowane do wyznaczania liczb kryterialnych.. Mimo to, z uwagi na nierównomierności i nieregularności powierzchni. kawern, konieczna jest weryfikacja empiryczna współczynników wymiany ciepła (Urbańczyk in., 2011) Dotychczasowe obliczenia geomechaniczne wykonywane w Polsce nie uwzględniały zjawisk termodynamicznych zachodzących w kawernie i ich wpływu na otaczający górotwór. Wyznaczając stan deformacyjno-naprężeniowy górotworu przyjmowano, że temperatura jest stała w czasie i zależy tylko od głębokości zgodnie z wartością gradientu geotermicznego. W świecie prace badawcze dotyczące wpływu zjawisk termodynamicznych na przemieszczenia, odkształcenia i naprężenia górotworu prowadzone są od kilku lat np. (Rokahr, 2008). Podjęto je nie tylko w związku z magazynowaniem gazu ziemnego lecz również magazynowaniem energii w postaci sprężonego powietrza, określanego w literaturze jako technologia CAES (Compressed Air Energy Storage) (Crotogino i in., 2001), (Crotogino.,. 2006).. W. pierwszej. kolejności. odpowiednie. obliczenia. numeryczne. wykonywano przy wykorzystaniu dwóch niezależnych symulatorów: pierwszy analizujący problemy termodynamiczne oraz drugi opisujący zagadnienia geomechaniczne. W pracach (PB Energy Storage Services, 2011), (Karimi–Jafari, i in., 2011), (Leuger i in, 2012) wykorzystano opracowany w USA przez PB Energy Storage Services i RESPEC Inc symulator SCTS (Salt Cavern Thermal Simulator) przeznaczony do symulacji zjawisk termodynamicznych i przepływu ciepła powstałych w trakcie eksploatacji kawern magazynowych wykonanych w złożach soli kamiennej (Nieland, 2004). Przykładowo w czasie symulacji komory o objętości ok. 155 tyś. m3, eksploatowanej w przedziale ciśnień od 5,5 do 10,3 MPa uzyskano zmiany wartości temperatury gazu około 30ºC. Obliczenia geomechaniczne wykonywano wykorzystując stosowane wcześniej komercyjne programy takie jak ABAQUS czy FLAC. W programach tych zmienność temperatury wprowadza się poprzez warunki brzegowe pierwszego i drugiego rodzaju tzn. jako wartości temperatury lub wielkości strumienia cieplnego na brzegu kawerny. Programem który jest dostosowany do rozwiązywania obu rodzajów zagadnień jest program LOCAS opracowany przez Brouard Consulting (Brouard i in. 2006). Program ten dedykowany jest do kawern solnych i uwzględnia następujące zjawiska: krótko i długotrwałe pełzanie soli, ogrzewanie solanki i rozchodzenie się ciepła, mikro-przesiąkanie solanki w oparciu o prawo Darcy, rozpuszczanie i rekrystalizację soli, adiabatyczne sprężanie i 8.

(10) rozprężanie oraz zniszczenie soli (zwiększenie przepuszczalności) wg wybranego spośród trzech kryterium. Obliczenia wykonane tym programem przedstawiono m.in. w pracach (Karimi–Jafari, i in., 2011) (Brouard i in. 2012). Program LOCAS był reklamowany na sympozjum SMRI i do niektórych wymienianych jego możliwości podejść należy sceptycznie. W pracy (Karimi–Jafari, i in., 2011) przedstawiono obydwa opisane powyżej podejścia. W pierwszym rozkład temperatur na ociosie kawerny uzyskany programem SCTS był wczytany do programu ABAQUS, w drugim zastosowano sprzężone termo-mechaniczne symulacje programem LOCAS. Obliczenia wykonano dla cylindrycznej kawerny (wys. 400 m, średnica 70 m) o objętości 1 482 000 m3 zakończonej kopułami, zlokalizowanej na głębokości 1400 m p.p.t. Symulowano dwa 35 letnie scenariusze eksploatacji: tradycyjny sezonowy oraz łączony sezonowy i wielocykliczny. W trakcie obliczeń nie uwzględniono wychłodzenia górotworu na skutek ługowania a temperatura gazu na końcu pierwszego napełniania była równa geotermalnej. W przypadku obliczeń programem SCTS symulowano sferyczny kształt komory a wartości ciśnień regulowane były wydajnością poboru i zatłaczania gazu. Uwzględniono również wychłodzenie górotworu na skutek ługowania. W podsumowaniu autorzy stwierdzają, że bardzo częste cykle napełniania i opróżniania komory nie pogarszają stateczności komory w porównaniu z jej eksploatację sezonową. Analiza rezultatów nie obejmowała jednak konwergencji i prędkości konwergencji kawerny.. 1.4. Cel i zakres pracy Niniejsza. praca. realizowana. jest. w. ramach. projektu. badawczego. pt.. „Termomechaniczna analiza procesu eksploatacji kawernowego podziemnego magazynu gazu” realizowanego na Wydziale Wiertnictwa Nafty i Gazu AGH, którego celem jest opracowanie. modelu. numerycznego. symulującego. stan. deformacyjno-naprężeniowy. górotworu z uwzględnieniem procesów termodynamicznych zachodzących w kawernie solnej. Model ten ma łączyć algorytmy dwóch wcześniej opracowanych programów numerycznych: • programu Kaga opracowanego w OBR Chemkop (Kunstman i in., 1980) przeznaczonego do wyznania rozkładów temperatury w kawernie i otaczającym ją górotworze, którego algorytm zostanie przedstawiony w rozdziale 4,. 9.

(11) • programu GEOSOLK opracowanego w OBR Chemkop i IGSMiE PAN przeznaczonego do wyznaczania stanu deformacyjno-naprężeniowego górotworu w otoczeniu kawerny magazynowej jednak bez uwzględnienia zjawisk termicznych. Praca stanowi pierwszy etap projektu badawczego i jej celem było wykonanie obliczeń dla wyidealizowanej cylindrycznej kawerny wykorzystując oba powyższe modele bez znaczących ingerencji w ich algorytmy. Do przeprowadzenia obliczeń skorzystałam z opracowanego już wcześniej programu konwertującego wyników pomiędzy programem Kaga opartym o metodę różnic skończonych przy wykorzystaniu siatek prostokątnych a programem GEOSOLK opartym o metodę elementów skończonych w użyciem siatek nieregularnych. Zakres przeprowadzonych w pracy obliczeń i analiz był więc następujący: • wykonanie obliczeń programem Kaga dla kawern o średnicy 50 m, wysokości 200 m których środek zlokalizowany był na głębokości 800, 1100 i 1400 m p.p.t. eksploatowanych przez okres 5 lat zgodnie z założonym scenariuszem, • wprowadzenie. uzyskanych. rozkładów. temperatur. górotworu. do. programu. GEOSOLK i wyznaczenie wartości przemieszczeń i odkształceń górotworu przy 4 różnych wariantach prawa pełzania Nortona których współczynniki określono na podstawie badań laboratoryjnych różnych rodzajów soli, • porównanie otrzymanych rezultatów z wynikami uzyskiwanymi dotychczasową wersją programu GEOSOLK. W sumie analiza objęła 32 modele numeryczne. Zaprezentowana metodyka obliczeń jest jakościowo zbliżona do analiz wykonanych programami SCTS i ABAQUS, wykorzystuje jednak krajowe programy komputerowe wykorzystywane wcześniej przy ocenie funkcjonowania kawern magazynowych gazu w KPMG (Kawernowy Podziemny Magazyn Gazu) Mogilno i KPMG Kosakowo. Dodatkowo w pracy zamieszczono wyniki badań laboratoryjnych, których celem było sprawdzenie poprawności prawa pełzania Nortona w warunkach odpowiadających eksploatacji komory magazynowej i weryfikacja jego współczynników. Ostatecznym celem pracy była wszechstronna analiza spadku pojemności kawerny magazynowej w czasie w zależności od parametrów termodynamicznych (temperatury wychłodzenia górotworu i zmian temperatury wywołanej sprężaniem i rozprężaniem gazu) oraz parametrów dodatkowych takich jak głębokość posadowienia kawerny, prędkość pełzania soli i wartość współczynnika określającego wpływ temperatury na prędkość tego procesu. Podstawową analizowaną wielkością była konwergencja. 10.

(12) względna kawerny czyli spadek jej objętości w czasie w stosunku do objętości początkowej. W pracy omówiono również rozkład przemieszczeń górotworu w otoczeniu kawerny. Konwergencja występuje w przypadku wszystkich wyrobisk w złożach soli, przy czym można ją rozpatrywać jako zjawisko pozytywne lub negatywne (Maj, 2009). Do pierwszego przypadku zaliczyć należy składowanie odpadów niebezpiecznych, kiedy to konwergencja gwarantuje szczelne zamknięcie wyrobisk oraz magazynowanie ropy naftowej i paliw, w których największym problemem jest rozługowywanie kawerny w trakcie wypychania surowca wodą lub solanką. W przypadku magazynowania gazu konwergencję należy jednoznacznie traktować jako zjawisko negatywne.. 11.

(13) 2. Teoretyczne podstawy modeli termodynamicznych kawern magazynowych gazu W rozdziale 1.2 wymieniono 4 zjawiska termodynamiczne, które mają wpływ na temperaturę a w konsekwencji na stan naprężeń i odkształceń w otoczeniu komory i jej konwergencję. Dotyczą one zmian zachodzących w magazynowanym gazie, wymiany ciepła między górotworem a medium ługującym lub gazem oraz przewodnictwem cieplnym górotworu. W algorytmie programu KAGA rozwiązywany jest numerycznie układ 7 równań różniczkowych które przedstawiono schematycznie na rysunku 2.1.. Strzałki na rysunku wskazują którego z czterech elementów komory magazynowej (magazynowany gaz, otwór udostępniający, brzeg kawerny, górotworu otaczający) dotyczy równanie.. Rys. 2.1. Równania rozwiązywane w programie KAGA. 12.

(14) Poniżej podano komplet równań w sformułowaniu całkowym, użyty w programie KAGA zaczerpnięty z opracowania (Urbańczyk i in. 1991). W kolejnych podrozdziałach podano równania dotyczące kolejno: procesów zachodzących w magazynowanym gazie, wymiany ciepła pomiędzy gazem a górotworem oraz zmiany temperatury w górotworze solnym.. 2.1. Procesy zachodzące w magazynowanym gazie. 2.1.1. Równanie stanu magazynowanego gazu Jeśli własności gazu odbiegają od gazu doskonałego, tradycyjnie równanie stanu gazu doskonałego poprawia się wprowadzając do niego czynnik ściśliwości Z. pv. =. RT Z. (2.1). gdzie: p – ciśnienie R – uniwersalna stała gazowa T – temperatura v – objętość molowa Czynnik ściśliwości Z nie jest stałą, ale dość skomplikowaną funkcją stanu gazu. Dlatego równanie powyższe nie jest wystarczające, a potrzebne jest równanie z czynnikiem Z wyrażonym explicite jako funkcja zmiennych stanu. Głównym składnikiem gazu ziemnego jest metan (na ogół powyżej 95%), a warunki termiczne panujące w kawernie magazynowej są dość odległe od punktu krytycznego. Z tego powodu dopuszczalne jest opisanie zachowania się gazu przy pomocy równania stanu Redlicha-Kwonga. Równanie to ma postać następującą:. p. =. RT a − 0.5 v−b T v (v + b ). (2.2). gdzie:. a, b – stałe, zależne od składu gazu. 13.

(15) Dla wyznaczenia stałych a i b dla danego składu gazu, przyjmuje się następujące reguły mieszania: n. ∑ yi y j aij. a =. (2.3). i, j =1. n. b =. ∑ yi bi. (2.4). i =1. gdzie:. yi – molowy udział i-tego składnika w mieszaninie gazowej. aij =. (. ). 0.25 Ω ai +Ω a j RTc1ij.5 (v i +v j ) 0.291−0.04(ωi +ω j ) = Ω bi. bi. Tcij. =. R Tci. (2.5). (2.6). pci. Tci yi + Tc j y j yi + y j. (2.7). ωi – czynnik acentryczny i-tego składnika, wg odpowiednich tabeli pci – ciśnienie krytyczne i-tego składnika wg odpowiednich tabeli. T ci – temperatura krytyczna i-tego składnika wg odpowiednich tabeli Ω ai – stała, klasycznie równa 0.42748 , jednak można poprawić dokładność równania,. biorąc dla każdego składnika wartość „poprawioną” wg odpowiednich tabeli, Ω bi – stała, klasycznie równa 0.08664 , również można dla danego składnika brać. wartość „poprawioną” wg odpowiednich tabeli.. Czynnik ściśliwości Z dla równania Redlicha-Kwonga można przedstawić następująco:. Z. =. pv RT. =. 1 b 1− v. − T. 1.5. a R  b v 1 +   v. (2.8). postać ta wynika bezpośrednio z równania (2.2). W postaci (2.2) równania Redlicha-Kwonga występuje objętość molowa. W praktyce dogodniejsze jest stosowanie zamiast niej gęstości gazu, co wymaga przeliczenia stałych R, a, b na odpowiednie jednostki. 14.

(16) p. =. RT 1 ρ. −. −b. T. 0.5. 1 ρ. a  1 + b ρ   . (2.9). 2.1.2. Równanie ciągłości gazu Zmiennymi opisującymi stan gazu są temperatura, gęstość i ciśnienie związane wzajemnie równaniem stanu oraz prędkość przepływu. Wyznaczyć je można rozwiązując układ różniczkowo-całkowych równań ciągłości, ruchu i energii. Pierwszym z nich jest równanie ciągłości. Wyraża ono zasadę zachowania masy. Bilansuje zmianę masy w danej objętości z wypływem masy przez powierzchnię brzegową tej objętości następująco: ∂ ∂t. ∫∫∫ ρ dV. =. Ω. − ∫∫ ρ v ⋅ n dS. (2.10). Γ. gdzie: t – czas v – prędkość przepływu gazu. ρ – gęstość gazu Ω – dowolny element kawerny lub otworu, wypełniony gazem Γ – powierzchnia zewnętrzna elementu Ω n – wektor jednostkowy, normalny do powierzchni Γ, skierowany na zewnątrz. dV – element całkowania po objętości dS – element całkowania po powierzchni. 2.1.3. Równanie ruchu gazu Równanie to, wyrażające zasadę zachowania pędu, wiąże zmianę pędu gazu zawartego w danej objętości ze strumieniem pędu przepływającym przez powierzchnię brzegową oraz działaniem sił powierzchniowych i objętościowych (w tym przypadku grawitacyjnych) następująco:. 15.

(17) ∂ ∂t. ∫∫∫ ρ v dV. = − ∫∫ ρ vv ⋅ n dS +. Ω. ∫∫ σ$ ⋅ n dS. Γ. Γ. +. ∫∫∫ ρ g dV. (2.11). Ω. gdzie: g – wektor przyspieszenia ziemskiego |g|=g. σˆ – tensor naprężeń w gazie pozostałe oznaczenia jak w równaniu (2.10). Tensor naprężeń dogodnie jest rozdzielić na naprężenia normalne (związane z ciśnieniem) i naprężenia styczne (opory lepkie): σ$. =. $ − p $I + Λ. (2.12). gdzie: p – ciśnienie $ – tensor naprężeń stycznych Λ $I – tensor jednostkowy. Aby wyznaczyć rozkład pędu wewnątrz zbiornika należy podać warunki na brzegu. W przypadku kawerny gazowej można wyróżnić dwa typy brzegów: • brzeg przez który możliwa jest jedynie wymiana ciepła, • brzeg przez który następuje przepływ gazu wraz z zawartą w nim energią wszystkich rodzajów. Do pierwszego typu brzegu zalicza się powierzchnię kawerny oraz boczną powierzchnię kolumny rur w otworze. Strumień pędu jest na nim równy zeru, tensor naprężeń określa się na podstawie równania oporów przepływu. Do drugiego typu brzegu zalicza się wlot rury na głowicy i wylot rury eksploatacyjnej w szyi kawerny. Na podstawie wydajności przepływu określa się na tym brzegu odpowiednie strumienie brzegowe. Jeżeli brak przepływu, równanie (2.11) traci dwie pierwsze całki zależne od v i przechodzi w równanie równowagi pneumostatycznej:. ∫∫ p n dS Γ. =. ∫∫∫ ρ g dV. (2.13). Ω. biorąc pod uwagę, że przy braku przepływu nie pojawiają się w gazie naprężenia styczne.. 16.

(18) Dla gazu w kawernie jedynie równanie (2.13) jest stosowalne, gdyż z uwagi na średnicę kawerny, pozostałe wyrazy równania (2.11) są do zaniedbania. Należy zwrócić uwagę, że stratyfikacja gazu w kawernie możliwa jest tylko przejściowo, podczas zatłaczania gorącego gazu. Podczas poborów, jak i podczas postojów, dochodzi do mieszania się gazu zawartego w kawernie, z uwagi na wymianę ciepłą z górotworem, w którym temperatura rośnie z głębokością.. 2.1.4. Równanie energii gazu Równanie energii bilansuje zmianę energii w danej objętości z wypływem energii przez powierzchnię brzegową tej objętości.. Różnica w bilansie pochodzi od ciepła, które. wypłynęło przez powierzchnię brzegową oraz od pracy wykonanej przez siły zewnętrzne na gazie zawartym w tej objętości. W postaci globalnej równanie przedstawia się następująco:. ∂ ρ edV =− ∫∫ ρ ev ⋅ ndS + ∫∫ (σˆ ⋅ v ) ⋅ ndS + ∫∫ Φ ⋅ ndS + L ∂ t ∫∫∫ Ω Γ Γ Γ. (2.14):. gdzie: e – energia właściwa gazu (całkowita) v – prędkość przepływu gazu. ρ – gęstość gazu. σˆ – tensor naprężeń w gazie, obejmujący ciśnienie i naprężenia styczne (lepkie) Φ – strumień cieplny L – praca mechaniczna wykonana na gazie Należy zwrócić uwagę na fakt, że w równaniu (2.14) pracę wykonaną przez gaz rozdzielono na pracę naprężeń (siły powierzchniowe) -. ∫∫ (σ$ ⋅ v) ⋅ n dS. i pozostałą pracę sił. Γ. objętościowych - L. Całkowita energia właściwa składa się z trzech składników: energii wewnętrznej, energii kinetycznej i energii potencjalnej w polu sił grawitacyjnych, co można wyrazić się wzorem: e = uRT +. v2 − gh 2. (2.15). 17.

(19) gdzie: u – zredukowana energia wewnętrzna gazu R – stała gazowa T – temperatura gazu v – prędkość przepływu gazu g – przyspieszenie ziemskie h – głębokość od powierzchni terenu w dół Energię potencjalną można częściowo z równania (2.14) wyrugować, wykorzystując równanie ciągłości. Jeśli zaś za powierzchnię Ω przyjąć całą kawernę, lub odcinek rury, można również pracę sił powierzchniowych ograniczyć do pracy przeciw ciśnieniu, gdyż na bocznej powierzchni nie występuje normalna składowa prędkości, albo naprężenia styczne. Wtedy pod całką powierzchniową można z kolei połączyć ciśnienie z energią wewnętrzną w jedną funkcję termodynamiczną - entalpię. Równanie (2.14) podobnie jak równania (2.10) i (2.11) wymaga określenia warunków na brzegu rozpatrywanego obszaru, aby wyznaczyć warunki panujące wewnątrz niego. Jak była mowa, można wyróżnić dwa typy brzegów : • brzeg przez który możliwa jest jedynie wymiana ciepła • brzeg przez który następuje przepływ gazu wraz z zawartą w nim energią wszystkich rodzajów Na brzegu pierwszego typu wymianę ciepła wyraża się następującym związkiem:. ∫∫ Φ ⋅ n dS. =. Z. ∫∫ α (T. z. − T ) dS. (2.15). Z. gdzie: α – współczynnik wymiany cieplnej Tz – temperatura w ośrodku sąsiadującym T – temperatura gazu. Dla brzegu drugiego typu przepływy energii określa się na podstawie równań ruchu.. 18.

(20) 2.1.5. Równanie oporu przepływu przez kolumnę rur. Równanie (2.11) z podstawieniem (2.12) wymaga określenia tensora naprężeń stycznych.. Pojawiają się one tam, gdzie przepływ gazu osiąga dużą prędkość, czyli. praktycznie rzecz biorąc wyłącznie w kolumnie rur. Zagadnienie można sprowadzić do określenia oporów przepływu przez tę kolumnę. Z różnych form opisu, najpraktyczniejsza jest formuła Darcy'ego, formalnie poprawna dla przepływu ustalonego, stosowana w praktyce także dla przepływu ciągłego, bez skoków prędkości. Wiąże ona opór lepki (naprężenia styczne) ze średnią prędkością przepływu i geometrią:. ∫∫ Λ$ ⋅ n dS. =. −. Z. v π λ D ρ v 2 dl ∫ vL 8. (2.16). gdzie: ρ – gęstość solanki ˆ – tensor naprężeń stycznych w gazie Λ. λ – bezwymiarowy współczynnik oporu v – średnia w przekroju prędkość przepływu v = |v| D – średnica hydrauliczna rury. Z – powierzchnia boczna rozpatrywanego odcinka rury (do całki po tej powierzchni redukuje się całka po powierzchni zamkniętej, gdyż na przekrojach poprzecznych rury nie pojawiają się naprężenia styczne) L – długość wycinka przepływu, tzn. rozpatrywanego odcinka rury n – wektor jednostkowy, normalny do powierzchni Z, skierowany zewnątrz dS, dl – elementy całkowania po powierzchni i po długości rury.. 19.

(21) 2.2. Wymiana ciepła pomiędzy gazem a górotworem. Wymianę ciepła między gazem a górotworem opisuje równanie:. ∫∫ λ ∇T g. g. ⋅ ndS =∫∫ α (Tk − Tg )dS. Z. (2.31). Z. gdzie:. Tg – temperatura górotworu Tk – temperatura wewnątrz kawerny Z – dowolny fragment powierzchni brzegowej zbiornika α – współczynnik wymiany cieplnej Kawernę trudno traktować jako naczynie kuliste, czy walcowe, jej powierzchnia boczna jest znacznie większa, niż kuli o analogicznej objętości, czy regularnego walca. Kształty kawern w większości przypadków wykazują znaczne nieregularności, dlatego współczynnik. α wyznaczony teoretycznie należy traktować raczej jako ograniczenie od dołu, a właściwą wartość dobrać modelując przebieg rzeczywistej eksploatacji kawerny magazynowej, metodą prób i błędów, jako punkt wyjścia przyjmując wartość półtora – dwukrotnie większą, zależnie od nieregularności kawerny. Na szczęście, w wielu przypadkach praktycznych, niezbyt dokładna znajomość współczynnika wymiany cieplnej ma stosunkowo niewielki wpływ na temperaturę magazynowanego gazu.. Przy stosunkowo wysokich współczynnikach α,. powyżej 20 W/m2K, dopływ ciepła z górotworu zaczyna być limitowany głównie przez jego przewodnictwo cieplne. Inaczej przebiega wymiana ciepła w rząpiu kawerny, gdzie osiadają części nierozpuszczalne a przestrzeń pomiędzy nimi („pory”) pozostaje wypełniona solanką. Również nad osadami, po pierwszym napełnieniu kawerny gazem, pozostaje parę metrów solanki rezydualnej (Urbańczyk i in. 2011).. 20.

(22) 2.3. Zmiany temperatury w górotworze solnym – równanie przewodnictwa cieplnego Stan górotworu opisany jest równaniem przewodnictwa przyrównującym zmianę ciepła w danej objętości do bilansu przewodzenia ciepła przez powierzchnię boczną tej objętości.. Równanie to w postaci globalnej przedstawia się następująco:. ∂ ∂t. ∫∫∫ ρ. g. C g Tg dV. Ω. =. ∫∫ λ. g. ∇Tg ⋅ n dS. (2.47). Γ. gdzie :. ρg – gęstość górotworu Cg – ciepło właściwe górotworu Tg – temperatura t – czas λg – współczynnik przewodnictwa cieplnego w górotworze Ω – dowolny fragment górotworu Γ – powierzchnia fragmentu Ω. n – wektor jednostkowy, normalny do powierzchni Γ, skierowany na zewnątrz dV – element całkowania po objętości dS – element całkowania po powierzchni. Równanie (2.47) pozwala wyznaczyć rozkład temperatur w obszarze, na którego brzegu znana jest temperatura albo jej gradient. Rozpatrując górotwór wokół zbiornika gazu ziemnego, można wziąć tak duży obszar, że przemiany termodynamiczne w zbiorniku nie będą miały wpływu na warunki panujące na brzegu modelu. Można więc przyjąć, że na brzegu modelu temperatura nie zależy od czasu a jedynie od głębokości, zgodnie ze stopniem geotermicznym: Tg. = const ( h). (2.48). gdzie:. h - głębokość od powierzchni terenu. 21.

(23) 2.4. Algorytm programu KAGA Program KAGA przeznaczony jest do wyznania gęstości, ciśnienia i temperatury magazynowanego gazu oraz rozkładu temperatury w górotworze otaczającym cylindryczną kawernę. Schemat działania programu jest następujący (Urbańczyk K. i in. 1991): 1. Program wczytuje z pliku tekstowego lub z klawiatury dane dotyczące wielkości kawerny i głębokości jej posadowienia, wielkość kroków aproksymacyjnych dla górotworu, oraz profil przewodnictwa cieplnego skał, ich pojemność cieplną oraz współczynniki wymiany ciepła między gazem a górotworem w otworze i w kawernie. Na tej podstawie tworzy model górotworu i otworu prowadzącego do kawerny, stosując aproksymację różnicową. 2. Wczytuje z pliku tekstowego lub klawiatury dane dotyczące stanu początkowego modelu tj. rozkładu temperatury w górotworze, otworze i kawernie. 3. W razie potrzeby modyfikuje krok czasowy i parametry iteracyjne obliczeń. 4. Wczytuje z pliku tekstowego lub z klawiatury scenariusz eksploatacji, składający się z etapów, które dotyczą zatłaczania gazu, poboru gazu, lub przerwy między tymi operacjami. Każdy z etapów eksploatacji opisany jest niezależnym segmentem: a) Etap zatłaczania symulowany przez segment RUZA określają: czas trwania etapu, wydajność zatłaczania i temperatura zatłaczanego gazu. b) Etap poboru symulowany przez segment RUPO opisany jest przez czas trwania etapu i wydajność poboru. c) Etap przerwy symulowany przez segment PRZERWA określony jest przez czas trwania etapu. Każdy etap eksploatacji opisuje kilka pętli iteracyjnych, z których zewnętrzna dotyczy kroku czasowego eksploatacji (Urbańczyk i in. 2011).. 2.4.1. Opis segmentu zatłaczania gazu RUZA. Celem segmentu jest wyznaczenie w każdym kroku czasowym gęstości i temperatury gazu na wlocie do komory. Pętle wewnętrzne dotyczą iteracji w otworze wydobywczym, kawernie i górotworze. Uproszczony schemat blokowy przedstawiony jest na rysunku 2.2.. 22.

(24) Aproksymacja czasu i aproksymacja głębokości dla przepływu gazu przez rurę wydobywczą są sprzężone w uwikłanym schemacie Leleviera. Wielkości iterowane w każdym odcinku rury to gęstość, temperatura i prędkość przepływu gazu.. Aby rozpocząć proces iteracyjny trzeba określić oprócz temperatury. zatłaczanego gazu jego gęstość. Potrzebne jest do tego ciśnienie głowicowe, które szacuje się na podstawie warunków początkowych.. Po przeiterowaniu wszystkich segmentów rury,. otrzymuje się temperaturę i gęstość gazu bezpośrednio wpływającego do kawerny, skąd wynika ciśnienie w jej stropie. Równocześnie ciśnienie to wynika z gęstości gazu w kawernie i jego temperatury. Różnica tego ciśnienia w stropie wynikającego z dwóch metod jego wyznaczenia stanowi poprawkę, która musi być mała, aby przejść do iteracji temperatury w górotworze. Jeśli poprawka ta nie jest mała, ani nie jest mała suma poprawek temperatury, konieczne jest ponowne przeiterowanie segmentów rury, przy czym ciśnienie głowicowe zmienia się o określoną część poprawki ciśnienia (podrelaksacja). Równania w górotworze iteruje się metodą kolejnych nadrelaksacji, aż suma poprawek jest dostatecznie mała. Zanim nastąpi przejście do następnego kroku czasowego, powtarza się naprzemian iteracje gazu i iteracje górotworu, aż różnice staną się małe. Ponieważ temperatura górotworu jest wolniej zmienna niż temperatura gazu, można górotwór iterować z krokiem czasowym większym niż krok czasowy iteracji gazu. W związku z tym że równanie Redlicha-Kwonga jest nieliniowe, również równania gazu w każdym segmencie rury muszą być iterowane. Nie jest jednak celowe stosowanie tu wysokiej dokładności, gdyż w następnej pętli iteracyjnej otworu warunki ulegną zmianie. Praktyka pokazała, że najefektywniej jest wykonywać każdorazowo 10 iteracji równań gazu.. 23.

(25) Rys. 2.2. Schemat blokowy segmentu opisującego zatłaczanie gazu w programie KAGA. 24.

(26) 2.4.2. Opis segmentu PRZERWA Etap ten symulaowany jest w oparciu o taki sam schemat iteracyjny jak w segmencie poprzednim. Jego schemat blokowy przedstawiony jest na rysunku 2.3. Wielkości iterowane w każdym odcinku rury to gęstość i temperatura.. Prędkość. przepływu gazu jest bliska zeru, a w związku z tym jest do zaniedbania, ze względu na jej wpływ na temperaturę. Podczas iteracji równań w odcinkach rury i w kawernie, temperaturę wylicza się z równania energii (z uwzględnieniem wymiany ciepła z górotworem), ponieważ zmiany temperatury zmieniają ciśnienie, zmienia się czynnik ściśliwości Z, co nieznacznie zmienia gęstość gazu. Równania iteruje się zawsze 10-krotnie. Praktyka pokazała, że nie jest jednak celowe iterowanie do wysokiej dokładności, gdyż w następnej pętli iteracyjnej otworu warunki ulegną zmianie. Konieczność iteracji wynika z nieliniowości równania Redlich-Kwonga. Pętlę z iteracjami gazu w kawernie i segmentach rury wydobywczej powtarza się, aż suma poprawek temperatury będzie dostatecznie mała. Równania przewodnictwa w górotworze iteruje się metodą kolejnej nadrelaksacji, aż suma poprawek jest dostatecznie mała. Zanim nastąpi przejście do następnego kroku czasowego, powtarza się naprzemian iteracje gazu i iteracje górotworu, aż różnice staną się małe. Temperatura gazu podczas przerwy zmienia się powoli, w tempie podobnym do zmian temperatury górotworu.. W związku z tym stan gazu iteruje się takim samym krokiem. czasowym jak temperatury w górotworze, dłuższym, niż przy zatłaczaniu gazu lub jego poborze.. 25.

(27) Rys. 2.3. Schemat blokowy segmentu opisującego etap przerwy w programie KAGA. 26.

(28) 2.4.3. Opis segmentu poboru gazu z kawerny RUPO. Etap ten symulaowany jest w oparciu o taki sam schemat iteracyjny jak schematy poprzednie. Schemat blokowy segmentu przedstawiony jest na rysunku 2.4. Wielkości iterowane w każdym odcinku rury to jego gęstość, temperatura i prędkość przepływu.. Aby rozpocząć proces iteracyjny trzeba określić oprócz temperatury gazu. wypływającego z kawerny jego gęstość i prędkość przepływu na wejściu z kawerny do rury wydobywczej. Temperatura i gęstość gazu wynikają z warunków w kawernie, wydajność szacuje się na podstawie głowicowej wydajności poboru. Po przeiterowaniu wszystkich segmentów rury, otrzymuje się temperaturę, gęstość gazu i jego szybkość przepływu przy głowicy.. Stąd wynika wydajność wypływającego gazu, którą porównuje się z zadaną. wydajnością poboru. Różnica tych wydajności stanowi poprawkę, która musi być mała, aby przejść do iteracji temperatury w górotworze. Jeśli poprawka ta nie jest mała, ani nie jest mała suma poprawek temperatury, konieczne jest ponowne przeiterowanie segmentów rury, przy czym szybkość przepływu na dole rury zmienia się o określoną część poprawki szybkości na głowicy (podrelaksacja). Równania przewodnictwa w górotworze iteruje się metodą kolejnej nadrelaksacji, aż suma poprawek będzie dostatecznie mała. Zanim nastąpi przejście do następnego kroku czasowego, powtarza się naprzemian iteracje gazu i iteracje górotworu, aż suma różnic stanie się dostatecznie mała tak jak dla etapu zatłaczania. Ponieważ temperatura górotworu jest wolniej zmienna niż temperatura gazu, można górotwór iterować z krokiem czasowym większym niż krok czasowy iteracji gazu.. 27.

(29) Rys. 2.4. Schemat blokowy segmentu opisującego etap poboru w programie KAGA. 28.

(30) 3. Teoretyczne podstawy modeli geomechanicznych kawern magazynowych gazu. 3.1. Równanie konstytutywne soli kamiennej Sól kamienna jest ośrodkiem wykazującym właściwości reologiczne. Charakterystyczną cechą takich ośrodków jest istotny wpływ czasu na wielkość odkształceń. Rysunek 3.1 przedstawia zależność pomiędzy odkształceniami i naprężeniami dla trzech procesów ściskania próbki przebiegających w czasie t1, t2, i t3, dla którego wraz ze wzrostem czasu. obciążenia rośnie wartość odkształceń uzyskiwanych przy danym poziomie naprężeń.. Rys. 3.1. Charakterystyki naprężeniowo-odkształceniowe ośrodka reologicznego (Ślizowski, 2006). W przypadku soli odkształcenia zależą nie tylko od naprężeń i czasu, lecz dodatkowo, w istotny sposób, od szeregu czynników takich jak temperatura, wilgotność czy budowa petrograficzna, a w szczególności wielkość kryształów (ziaren). Właściwości reologiczne powodują, że nawet w przypadku stałych obciążeń odkształcenie przyrasta w czasie, co przedstawia rysunek 3.2 na którym można wyróżnić można następujące fazy odkształceń:. • odkształcenie początkowe (sprężysto-plastyczne) powstające bardzo szybko dlatego zakłada się, że następuje natychmiastowo;. • pełzanie pierwotne w której prędkość pełzania maleje; • pełzanie stacjonarne tzn. w którym prędkość przyrostu odkształceń jest stała; • faza końcowa, w której przy dużych naprężeniach przekraczających 50% wytrzymałości doraźnej, pełzanie ulega gwałtownemu przyspieszeniu prowadzącemu. 29.

(31) do niszczenia materiału, przy mniejszych zaś dalej pozostaje stałe a poniżej 20 % zanika (Ślizowski, 2006).. Rys. 3.2. Fazy odkształceń soli kamiennej. W przypadku zmiany obciążeń przebieg procesu pełzania zależy od tego czy obciążenie maleje czy rośnie, co obrazuje rysunek 3.3. Jeśli w czasie t1 obciążenie zostaje zwiększone, pojawia się efekt pełzania pierwotnego a później prędkość pełzania spada do wartości pełzania stacjonarnego odpowiadającego nowym wartościom naprężeń. Jeśli natomiast w czasie t1 zmniejszymy obciążenie pełzanie początkowo zanika a potem stopniowo przyśpiesza do nowej prędkości pełzania stacjonarnego.. Rys. 3.3. Przebieg pełzania w przypadku zmiany obciążeń. 30.

(32) W związku z powyższym najwięcej uwagi poświęca się pełzaniu stacjonarnemu. Formułując równanie konstytutywne soli kamiennej (związek pomiędzy naprężeniami i odkształceniami) w większości modeli przyjmuje się że całkowite odkształcenie jest sumą odkształceń sprężystych i pełzania:. dε = dε el + dε v. (3.1). gdzie: dε - całkowity przyrost odkształceń dε el - przyrost odkształceń sprężystych (elastycznych) dε v - przyrost odkształceń lepkich (pełzania). W modelach przyjmuje się dodatkowo że sól kamienna jest izotropowa a zmiana stanu odkształcenia jest funkcją tylko aktualnego stanu zmiennych niezależnych tzn. napr ężenia, odkształcenia, czasu i temperatury (Ślizowski, 2006). Założenie superpozycji odkształceń umożliwia osobny opis każdego rodzaju odkształceń.. 3.1.1. Odkształcenia sprężyste. Opis odkształceń sprężystych soli oparty jest najczęściej na liniowym prawie Hooke’a i założeniu izotropii materiału. Prowadzi to do znanych zależności (Ślizowski, 1988):. ε ijel =. 1+ν 3ν σ ij − σ m δ ij E E. (3.2). 1 σm 3K. (3.3). εm =. gdzie:. ε ijel - składowe tensora odkształceń, σ ij - składowe tensora naprężeń,. 31.

(33) sij = σ ij − σ m ⋅ δ ij. eij = ε ij − ε m ⋅ δ ij. σ m = 13 (σ 11 + σ 22 + σ 33 ) - naprężenie średnie ε m = 13 (ε 11 + ε 22 + ε 33 ) - odkształcenie średnie 1 dla i = j. δ ij =. 0 dla i ≠ j. - delta Kroneckera. E – moduła Younga (odkształcalności podłużnej). ν - współczynnik Poissona G=. E - moduł Kirchhoffa (odkształcalności postaciowej) 2(1 + ν ). K=. E - moduł ściśliwości 3(1 − 2ν ). W niniejszej pracy przyjęto wartości E=8 GPa i ν=0,2.. 3.1.2. Odkształcenia pełzania Prędkość składowych tensora odkształceń pełzania wyraża zależność:. dε ijv dt. = 32 ⋅. dε efv dt. ⋅. sij. (3.4). σ ef. gdzie:. σ ef =. 2 2. (σ 11 − σ 22 ) 2 + (σ 22 − σ 33 ) 2 + (σ 33 − σ 11 ) 2 + 6(σ 12 + σ 23 + σ 31 ). ε ef =. 2 3. (ε 11 − ε 22 ) 2 + (ε 22 − ε 33 ) 2 + (ε 33 − ε 11 ) 2 +. 2. 3 2. 2. 2. 2. 2. 2. (ε 12 + ε 23 + ε 31 ). Niewiadomą w powyższym związku jest prędkość odkształceń efektywnych pełzania. dε efν / dt . Prawem pełzania nazywamy formułę opisującą prędkość odkształceń efektywnych pełzania jako funkcję naprężeń efektywnych, odkształceń efektywnych, czasu i temperatury oraz ewentualnie takich czynników jak parametry strukturalne, wilgotność i naprężenie średnie (Ślizowski, 2006):. 32.

(34) dε efν dt. = f (σ ef , ε νef , t , T , S ). (3.5). Stosowane są różne postacie funkcji f (Lindner, Brady 1983), (Lux, Hausermann 1983), (Lux, Rokahr 1984), (Langer 1984), (Hunshe, Hampel A., 1999), (Lux, Düsterloh, Hou 2001), (Flac 1998), (Maxwell 1984). W niniejszej pracy opis pełzania ograniczono do pełzania stacjonarnego przyjmując najczęściej stosowane prawo pełzania Nortona: dε ef dt. = Ae. −. Q RT. σ ef n. (3.6). Q – energia aktywacji, R – 8,3144 J⋅mol-1K-1 – stała gazowa, T – temperatura [K], A, n – stałe współczynniki. Wpływ. parametrów. strukturalnych. uwzględniono. rozpatrując. różne. wartości. współczynników A i Q/R przedstawione bliżej w rozdziale 5.3.. 3.2. Temperatura i naprężenie pierwotne w górotworze W świetle pomiarów in-situ temperatura pierwotna górotworu solnego zmienia się liniowo wraz z głębokością według wzoru (Ślizowski, Serbin, 2012):. T =α ⋅H + β. (3.7). gdzie:. α – gradient temperatury β – stały współczynnik H – głębokość [m p.p.t] Wartość gradientu temperatury jest różna w zależności od rozpatrywanej struktury solnej co przedstawia tabela 3.1.. 33.

(35) Tabela 3.1. Gradient temperatury α i współczynnik β stosowane do wyznaczenia temperatury pierwotnej górotworu Region Środkowopolskie wysady solne Wyniesienie Łeby Monoklina przedsudecka. Gradient temp. α [K/m] Współczynnik β 0,03-0,034 285 0,01 283 0,027 285. Najwyższą wartość przekraczającą 0,03 K/m gradient temperatury osiąga w wysadach środkowo-polskich. Nieco niższa jest wartość na monoklinie przedsudeckiej 0,027 K/m,. natomiast zdecydowanie najniższa jest wartość w pokładzie solnym na Wyniesieniu Łeby tj. 0,01 K/m (Ślizowski, Urbańczyk, 2011). Konsekwencją tego są znaczne różnice temperatury górotworu, która na głębokości 1000 m w przypadku wysadów wynosi 43°C a na Wyniesieniu Łeby 20°C. Powoduje to znaczne różnice pełzania górotworu gdyż jak wynika z analiz prawa pełzania w przypadku temp. 20°C prędkość pełzania jest około czterokrotnie niższa niż przy temp. 40°C. Problem naprężeń pierwotnych w soli w odróżnieniu od innych skał jest najmniej dyskusyjny. Powszechnie przyjmuje się, że wszystkie składowe naprężeń pierwotnych są sobie równe i maja wartość ciśnienia litostatycznego równego ciężarowi warstw nadległych. Założenie to jest przyjmowane ze względu na silne właściwości reologiczne górotworu solnego.. 3.3. Algorytm programu GEOSOLK Program GEOSOLK przeznaczony do rozwiązywania zagadnień symetrii osiowej ośrodka ciągłego zgodnie z teorią małych przemieszczeń i odkształceń metodą elementów skończonych. Rozwiązanie zagadnień sprężysto-plastyczno-lepkich realizowane jest poprzez iteracyjne poprawianie rozwiązania sprężystego. W przypadku sprężystym proces obliczeń realizowany jest w 3 etapach: 1. generacja siatki elementów skończonych, 2. złożenie macierzy sztywności struktury i jej diagonalizacja metodą eliminacji Gaussa z zapamiętaniem kolejności przeprowadzonych obliczeń,. 34.

(36) 3. wygenerowanie wektora sił węzłowych pochodzących od naprężeń i odkształceń początkowych, sił masowych i obciążeń zewnętrznych modelu oraz wykonanie na nim kolejnych operacji zapamiętanych w pkt. b i wyznaczenie rozwiązania. Zastosowany w systemie algorytm wykorzystujący pasmowość macierzy sztywności i rugujący stopnie swobody pojawiające się po raz ostatni, znacznie zmniejsza pojemność pamięci operacyjnej, niezbędnej do wykonania obliczań Do. uwzględnienia. nieliniowości. fizycznej. (odkształceń. plastyczno-lepkich). wykorzystano metodę odkszłceń początkowych w której odkształcenia plastyczno-lepkie traktowane są jako pozorne odkształcenia początkowe. Celem procesu iteracyjnego jest wyznaczenie poprawki, przy której nowy stan naprężeń jest taki, że kolejna poprawka jest jej dostatecznie bliska. Inaczej mówiąc, chodzi o doprowadzenie do sytuacji, w której wartość aktualnie obliczanej poprawki jest, z zadaną dokładnością, równa poprzednio podstawionej (Ślizowski J. 2006). Program GEOSOLK realizuje obliczenia zgodnie ze schematem przedstawionym na rysunku 3.4. Proces obliczeniowy rozpoczyna się od obliczeń wstępnych obejmujących utworzenie siatki elementów skończonych, wyznaczenie macierzy sztywności oraz obciążeń wynikających naprężeń pierwotnych górotworu. Zasadniczą część obliczeń programem GEOSOLK obejmują zaznaczone na schemacie 3.4 pętle: 1. Pętla zewnętrzna po etapach opisujących scenariusz eksploatacji komory tj. napełniania, postój i opróżnianie. 2. Pętla wewnętrzna dotycząca podetapów na które dzielimy dany etap. 3. Pętla iteracyjna w której obliczane są poprawki od nieliniowości równania konstytutywnego górotworu solnego tzn. uwzględniane jest jego pełzanie.. 35.

(37) Rys. 3.4. Schemat blokowy programu GEOSOLK. 36.

(38) W aktualnej wersji modelu oprócz obciążeń zewnętrznych modelu wyliczonych na każdym podetapie do każdego elementy wczytywana jest również wartość temperatury obliczona programem KAGA. Ze względu na to, że program KAGA i GEOSOLK mają różne siatki opracowano program interpolujący wartości temperatur w środkach elementów skończonych na podstawie wartości temperatur uzyskanych w najbliższych węzłach siatki programu KAGA. W pierwszej kolejności program znajduje prostokątne oczko programu KAGA w którym leży dany środek elementu skończonego siatki GEOSOLK. Następnie wartość temperatury w środku elementu o współrzędnych (x,y) interpolowana jest z wartości temperatur z wierzchołków oczka według schematu zamieszczonego na rysunku 3.5..  x − x   y2 − y   ⋅   T ( x, y ) = T ( x1 − y1 ) ⋅  2  x 2 − x1   y 2 − y1   x − x1   y 2 − y   ⋅   + T ( x 2 , y1 ) ⋅   x 2 − x1   y 2 − y1   x − x   y − y1   ⋅   + T ( x1 , y 2 ) ⋅  2  x 2 − x1   y 2 − y1   x − x1   y − y1   ⋅   + T ( x 2 , y 2 ) ⋅   x 2 − x1   y 2 − y1 . Rys. 3.5. Schemat obliczenia temperatury dla środków elementów generowanych w programie GEOSOLK. 37.

(39) 4. Laboratoryjne badania prędkości pełzania soli Celem badań było sprawdzenie funkcjonowania prawa pełzania Nortona w warunkach wieloetapowych procesów pełzania przy zmiennych temperaturach i obciążeniach. Badaniami objęto 3 rodzaje soli kamiennej drobno, średnio i grubokrystalicznej przedstawione odpowiednio na fotografiach 4.1a), 4.1b), 4.1c).. a). b). c). Rys. 4.1. Próby soli kamiennej wykorzystane do badań laboratoryjnych a) drobno-krystalicznej, b) średnio-krystalicznej, c) grubo-krystalicznej. W tabeli 4.1. zamieszczono scenariusze zmian obciążenia i temperatury odpowiednich prób soli w trakcie badań laboratoryjnych. Rysunki 4.2, 4.3 i 4.4 przedstawiają krzywe pełzania oraz prędkość pełzania soli odpowiednio drobno, średnio i grubokrystalicznej dla sześciu etapów zmian warunków ciśnienia i temperatury przedstawionych w tabeli 4.1.. 38.

(40) Tabela 4.1. Scenariusz zmian wartości ciśnienia i temperatury w trakcie badań laboratoryjnych prób soli kamiennej Rodzaj soli. Drobnokrystaliczna. Średniokrystaliczna. Grubokrystaliczna. P=10 [MPa] T=30°. 10. Nr etapu. Ciśnienie [MPa]. Temp. [ºC]. 1 2 3 4 5 6 1 2 3 4 5 6 1 2 3 4 5 6. 0 10 12,5 10 10 12,5 0 12,5 10 10 12,5 7,5 0 7,5 7,5 10 7,5 12,5. 30 30 30 50 30 60 30 30 50 30 30 60 30 30 50 30 30 60. P=12,5 [MPa] T=30°. P=10 [MPa] T=50°. P=10 [MPa] T=30°. P=12,5 [MPa] T=60°. 0.1 0.09. 9 0.08 8. 0.07. 0.05 6. ε [‰]. 0.04 0.03. 5. 0.02 4 0.01 3. szybkość ε [‰/dobę]. 0.06. 7. 0 -0.01. 2. -0.02 1 -0.03 0. -0.04 0. 10 20 30 40 50 60 70 80 90 100 110 120 130 140 150 160 170 180 190 200 210 220 230 240. czas [doby]. Rys. 4.2. Krzywe pełzania drobno-krystalicznej soli kamiennej. 39.

(41) P=12,5 [MPa] T=30°. 60. P=10 [MPa] P=10 [MPa] T=50° T=30°. P=12,5 [MPa] T=30°. P=7,5 [MPa] T=60°. 1.6 1.5. 56. 1.4 1.3. 52. 1.2. 48. 1.1 1. 44. 0.8 0.7. 36. 0.5. 28. 0.4 0.3. 24. 0.2 0.1. 20. 0 16. -0.1 -0.2. 12. -0.3. 8. -0.4 -0.5. 4. -0.6. 0. -0.7 0. 10 20 30 40 50 60 70 80 90 100 110 120 130 140 150 160 170 180 190 200 210 220 230 240. czas [doby]. Rys. 4.3. Krzywe pełzania średnio-krystalicznej soli kamiennej. P=7,5 [MPA] T=30°. 8. P=7,5 [MPA] T=50°. P=10 [MPA] T=30°. P=7,5 [MPA] P=12,5 [MPA] 0.1 T=30° T=60° 0.09. 7 0.08 0.07. 6. 0.05 0.04. 4. 0.03 3. szybkość ε [‰/dobę]. 0.06 5. ε [‰]. ε [‰]. 0.6 32. szybkość ε [‰/dobę]. 0.9 40. 0.02 0.01. 2. 0 1 -0.01 -0.02. 0 0. 10 20 30 40 50 60 70 80 90 100 110 120 130 140 150 160 170 180 190 200 210 220 230 240. czas [doby]. Rys. 4.4. Krzywe pełzania grubo-krystalicznej soli kamiennej. 40.

(42) Prędkości pełzania poszczególnych prób wyznaczono dla przedziałów czasowych w których proces pełzania był stacjonarny. Wykorzystano dwie metody: metoda 1 jako średnią wartość prędkości pełzania oznaczonych na rysunkach linią czerwoną oraz metoda 2 jako nachylenie stycznej do krzywej pełzania oznaczonej na rysunkach linią niebieską. Wyznaczone prędkości pełzania i przedziały aproksymacji zestawiono w tabeli 4.2.. Tabela 4.2. Prędkość pełzania stacjonarnego uzyskana przy naprężeniach i temperaturach dla których wyznaczono wartości współczynników prawa pełzania Nortona. Nr próbki. Ciśnienie [MPa]. Temp. [K]. Przedział aproksymacji [dni]. 1. Sieroszowice (sól średniokrystaliczna) 2. Mogilno (sól grubokrystaliczna) 3. Kosakowo (sól drobnokrystaliczna). 10 10 12,5 7,5 10 7,5 12,5 10 10. 323 303 303 323 303 303 303 323 303. 105-118 136-154 176-199 94-115 147-170 178-199 101-118 127-154 176-199. Prędkość pełzania stacjonarnego [‰/dobę] Metoda 1 (przeciętna prędkość). Metoda 2 (nachylenie stycznej). 0,19475 0,06256 0,18047 0,00800 0,00388 0,00290 0,00648 0,00618 0,00186. 0,18466 0,06125 0,18303 0,00728 0,00408 0,00276 0,00581 0,00565 0,00198. Różnica względna. 5,2 2,1 1,4 9,0 5,3 4,8 10,3 8,6 6,7. Wyznaczając wartości współczynników prawa pełzania skorzystano z metody zaprezentowanej w pracy (Ślizowski K. i in., 2005). Metoda to polega na analizowaniu poszczególnych próbek podczas zmienno-naprężeniowych i zmienno-temperaturowych testów pełzania. Taką metodykę przyjęto planując przebieg eksperymentów. Metoda pozwala na zupełne wyeliminowanie wpływu budowy petrograficznej poszczególnych próbek, wymaga jednak dłuższego czasu trwania poszczególnych testów. Wartości współczynników Q/R i n wyrażają się wzorami:. dε efvs (1) Q / R = ln. dt. ε. vs ef ( 2 ).  1 1 /  −   T2 T1 . (4.1). dt. 41.

(43) gdzie:. ε efvs (1). - prędkość odkształceń osiowych w temperaturze T1,. dt. ε efvs ( 2 ). - prędkość odkształceń osiowych w temperaturze T2.. dt. ε efvs(1) σ ef (1) n = ln vsdt / ln. σ ef ( 2 ). ε ef ( 2 ). (4.2). dt. gdzie:. ε efvs (1) dt. ε efvs ( 2 ) dt. - prędkość odkształceń osiowych przy ciśnieniu σ ef (1) - prędkość odkształceń osiowych przy ciśnieniu σ ef ( 2). Wartość współczynnika A obliczono zgodnie z wzorem:. ε efvs ( j ) dt. Aj =. −. e. Q RT( j ). (4.3). ⋅ σ efn ( j ). gdzie: dε j dt. , T j ,σ j. są odpowiednio prędkością pełzania, temperaturą i naprężeniem. efektywnym próbek na konkretnym etapie eksperymentu. Wartości współczynników prawa pełzania A, Q/R i n wyznaczone zgodnie z wzorami (4.1), (4.2), (4.3) dla soli drobno, średnio i grubo-krystalicznych przedstawia tabela 4.3. Przeprowadzone eksperymenty laboratoryjne wskazują że wielkość kryształów soli nie wpływa się na prędkość jej pełzania i przykładowo sole drobno-krystaliczne zachowują się podobnie do soli kryształowych. Badane sole średnio-krystaliczne wykazały się natomiast znaczenie wyższą prędkością pełzania, jednak z całą pewnością uzyskanych rezultatów nie można uogólniać. Generalnie można stwierdzić, że sól kamienna charakteryzuje się kilku lub nawet kilkunastokrotnym. 42.

(44) zróżnicowaniem prędkości pełzania i w związku z powyższym w dalszej pracy postanowiono rozważać dwie skrajne prędkości pełzania dopasowane do soli szybko i wolno pełznących.. Tabela 4.3. Współczynniki prawa pełzania Nortona uzyskane przy zmianach naprężeń i temperatury Nr próbki 1. 2. 3. Współczynnik prawa Nortona A Q/R n A Q/R n A Q/R n. Metoda 1. Metoda 2. 103,09 5556,94 4,75 4911,85 4963,16 1,01 0,18 3888,48 3,77. 41,82 5400,02 4,91 1124,21 4745,94 1,36 0,507 4840,95 4,56. Uzyskane podczas badań wartości współczynnika temperaturowego Q/R nie przekraczają 5000 K (patrz. tabela 4.3) tj. są niższe od wartości znanych z literatury (Hunshe, 1993), (Ślizowski, 2006). Z tego względu rozważano dalej dwie wartości współczynnika Q/R równe 5000K i 6500K. Również niższe niż w literaturze uzyskano wartości współczynnika n w związku z tym w dalszej części rozpatrzono dwie wartości tj. n=4 i n=5. Od strony jakościowej należy stwierdzić, że prawo Nortona poprawnie opisuje właściwości soli kamiennej zakresie zjawiska pełzania. Świadczą o tym stosunkowo niewielkie różnice względne pomiędzy wartościami teoretycznymi i rzeczywistymi przedstawione w tabeli 4.4. dla etapów eksperymentów, które nie były brane pod uwagę do wyznaczania współczynników.. Tabela 4.4. Rzeczywiste i teoretyczne prędkości pełzania uzyskane dla naprężeń i temperatury w których nie wyznaczano współczynników prawa pełzania Metoda 1 Nr próbki. Ciśnienie [MPa]. Temp. [K]. 1 2 3. 12,5 7,5 10. 303 303 303. Metoda 2. Rzeczywista. Teoretyczna. Różnica względna. Rzeczywista. Teoretyczna. Różnica względna. 0,1672 0,0033 0,0019. 0,1805 0,0029 0,0028. 8,0 12,1 47,4. 0,1594 0,0032 0,0020. 0,1830 0,0028 0,0021. 14,8 12,5 5,0. 43.

(45) 5. Zakres obliczeń numerycznych. 5.1. Warunki brzegowe i dyskretyzacja obszaru dla programu KAGA Na rysunku 5.1. przedstawiono siatkę aproksymacyjną wykorzystaną do obliczeń programem KAGA. 0 100 200 300 400 500 600 700 800 900 1000 1100 1200 1300 1400 1500 1600 1700 1800 1900 2000 0. 50. 100. 150. 200. 250. Rys. 5.1. Model obliczeniowy programu KAGA. 44.

(46) Ze względu na możliwości programu KAGA kształt modelowanych komór magazynowych został uproszczony do kształtu walca. Średnice osiowo-symetrycznych komór wynosiły 50 m a wysokość 200 m, co odpowiada objętości niemal 392.7 tys. m3, czyli podobnej do objętości, jaką mają przeciętne komory KPMG Mogilno (Ślizowski, Urbańczyk, 2011). Obliczenia wykonano dla trzech głębokości posadowienia środka komory tj. 800, 1100 i 1400 m p.p.t. Model obliczeniowy w programie KAGA składał się z 3160 elementów. Zasięgi poziomy i pionowy dla wszystkich głębokości były jednakowe i wynosiły 200 m od osi symetrii i 2000 m od powierzchni terenu. Na brzegu modelu przyjmuje się stałą temperaturę, wynikającą z pierwotnego rozkładu temperatury. Na konturze kawerny przyjmuje się warunek trzeciego rodzaju – wymianę ciepła z kawerną proporcjonalną do różnicy temperatur między górotworem a gazem w kawernie. Podobnie przyjmuje się wymianę ciepła w otworze, pomiędzy gazem a górotworem. Temperatura górotworu na głębokości środków kawern: 800, 1100 i 1400 m p.p.t. wyliczona zgodnie z wzorem (3.13) wynosiła odpowiednio 38.3°C, 47.3°C, 56.3°C. Do obliczeń przyjęto następujące wartości współczynników: •. przewodnictwa cieplnego górotworu: 6,5 W/mK,. •. przewodnictwa temperaturowego (dyfuzywność) górotworu: 0,42×10-5 m2/s,. •. wymiany ciepła w kawernie: 58,15 W/m2K,. •. wymiany ciepła w otworze: 5,82 W/m2K.. Wartości współczynników przewodnictwa cieplnego i przewodnictwa temperaturowego wyznaczono w oparciu o badania laboratoryjne soli z wysadu Mogilno. Współczynniki wymiany ciepła pomiędzy górotworem a gazem oszacowano na podstawie kalibracji modelu KAGA w oparciu o pomiary in-situ. Znaczne różnice pomiędzy współczynnikami wynikają z nieregularności kształtu kawerny (pozorne zwiększenie powierzchni wymiany i obecności płynu izolującego pomiędzy kolumnami rur.. 45.

(47) 5.2. Warunki brzegowe i dyskretyzacja obszaru dla programu GEOSOLK Rysunek 5.2. przedstawia podział modelowanego obszaru na elementy skończone. 700. 800. 900. 1000. 1100. 1200. 1300. 1400. 1500 0. 50 100 150 200 250. Rys. 5.2. Model obliczeniowy programu GEOSOLK Kształt i rozmiary analizowanej komory był identyczny jak w programie KAGA. Siatka elementów skończonych wykorzystana w programie GEOSOLK zbudowana jest z 1360 węzłów. Zasięg poziomy modelu wynosił 250 m, tj. 200 m od brzegu komory, zasięgi pionowe wynosiły po 250 m od stropu i spągu komory.. 5.3. Rozpatrywane prawa pełzania Opierając się na wynikach badań laboratoryjnych przedstawionych w rozdziale poprzednim i studiach literaturowych analizę konwergencji kawerny postanowiono przeprowadzić dla 4 praw pełzania Nortona różniących się wartością współczynników. Rozważono sole charakteryzujące się małą i dużą wrażliwością na zmiany temperatury. W. 46.

(48) pierwszym przypadku przyjęto, że współczynnik Q/R wyniesie 5000K a w drugim 6500K. Jak łatwo sprawdzić przyjęcie niższego współczynnika oznacza że prędkość pełzania w temperaturze 50°C jest 4.9 razy większa niż w temp 20°C natomiast przy przyjęciu Q/R = 6500 K różnica wynosi aż 7.9. Dla obu wartości współczynnika Q/R rozpatrzono dwie wartości współczynnika potęgowego n równego 4 lub 5. W tym przedziale mieści się najprawdopodobniej rzeczywista wartość współczynnika dla górotworu solnego. Rozpatrzono również sole szybko i wolno pełznące poprzez dopasowanie współczynnika A. Dla soli wolno pełznących przyjęto taką jego wartość by w temp. 50°C i naprężeniu 10 MPa prędkość odkształceń wyniosła 0.004‰/dobę natomiast dla soli szybko pełznących 0.032‰/dobę. Zestaw wszystkich czterech kompletów współczynników prawa pełzania Nortona wraz symbolem modelu zawiera tabela 5.1.. Tabela 5.1. Wartości współczynników prawa pełzania Nortna przyjęte do obliczeń numerycznych Model. A [‰]. Q/R [K]. n. dε/dt* [‰/dobę]. 1a. 2.113. 5000. 4. 0.004. 1b. 1.69. 5000. 5. 0.032. 2a. 219.65. 6500. 4. 0.004. 2b. 175.72. 6500. 5. 0.032. * prędkość wyliczona dla naprężeń efektywnych σef=10 MPa i temperatury T=323 K Uzyskane krzywe pełzania przedstawiano na rysunkach 5.3, 5.4 i 5.5. Na rysunkach tych przedstawiono również znane z literatury prawo BGRa (będące formalnie prawem Nortona) i BGRb (superpozycja dwóch praw Nortona) (Hunsche 1993, Hunsche Hampel 1999) oraz prawo Nortona ze współczynnikami wyznaczonymi w Instytucie Gospodarki Surowcami Mineralnymi i Energią PAN (Ślizowski i in. 2011).. 47.

(49) 1. prędkość εef [‰/dobę]. 0.1. 0.01. 0.001. Prawa pełzania. 0.0001. BGRa BGRb 1a 1b 2a 2b IGSMiE PAN. Wyniki testów pełzania Sieroszowice Mogilno Kossakowo. 1E-005 4. 6. 8. 10. 12. 14. 16. 18. 20. σef [MPa]. Rys. 5.3. Krzywe pełzania soli uzyskane w temperaturze 30ºC 10. prędkość εef [‰/dobę]. 1. 0.1. 0.01. Prawa pełzania 0.001. BGRa BGRb 1a i 2a 1b i 2b IGSMiE PAN. Wyniki testów pełzania Sieroszowice Mogilno Kossakowo. 0.0001 4. 6. 8. 10. 12. 14. 16. 18. 20. σef [MPa]. Rys. 5.4. Krzywe pełzania soli uzyskane w temperaturze 50ºC. 48.

(50) 1. prędkość εef [‰/dobę]. 0.1. 0.01. Prawa pełzania 0.001. Wyniki testów pełzania Sieroszowice Mogilno Kossakowo. BGRa BGRb 1a 1b 2a 2b IGSMiE PAN. 0.0001 292 294 296 298 300 302 304 306 308 310 312 314 316 318 320 322 324. T [K]. Rys. 5.5. Krzywe pełzania soli uzyskane przy ciśnieniu 10 MPa. 5.4. Wychłodzenie górotworu podczas ługowania i scenariusz eksploatacji kawerny magazynowej. Czas procesu ługowania kawern przeznaczonych na magazyny gazu może wynosić nawet kilka lat i uzależniony jest od uzyskiwanych stężeń. Generalnie jest jednak krótszy od niż w przypadku kawern eksploatacyjnych, gdzie pozyskiwana solanka musi osiągnąć stężenie nasycenia. W przypadku KPMG Mogilno uzyskiwana nienasycona solanka przepompowywana była do kopalni Mogilno zlokalizowanej na tym samym wysadzie. W analizowanej w pracy sytuacji nie rozpatrywano fazy ługowania komory, tzn. przyrostu w czasie jej wymiarów. W opracowanym modelu numerycznym przyjęto ostateczny kształt kawerny, dla którego symulowano konwergencję dla 4 lat postoju kawerny wypełnionej solanką co będzie widoczne na rysunkach przedstawionych w rozdziale 7.7.1. Różnice pomiędzy rozpatrywanymi dalej wariantami polegają na przyjęciu różnych pól temperatur początkowych górotworu w jej otoczeniu. Wielkość wychłodzenia górotworu ma bowiem oprócz efektów termodynamicznych zachodzących w kawernie istotny wpływ na 49.