Analiza cieplna procesu ochładzania górotworu przy uwzględnieniu ruchu wód podziemnych

(1)

ZESZYTY N A U K O W E P O L IT E C H N IK I Ś L Ą S K IE J Seria: E N E R G E T Y K A z. 124

1995 N r kol. 1278

Adam F IC , Ja n S K Ł A D Z IE Ń

ANALIZA CIEPLNA PR O C ESU OCHŁADZANIA

GÓROTWORU PRZY UW ZGLĘDNIENIU RUCHU WÓD PODZIEMNYCH

S tr e sz c z e n ie . W pra cy przedstaw iono sposób a n a lizy cieplnej proce

su ochładzania i wstępnego zam ra ża n ia górotw oru z uw zględnieniem ru c h u w ody g ru n to w e j o znanym ro zkła d zie prędkości oraz metodę obliczeń pola prędkości te j wody. Zaprezentow ano re z u lta ty obliczeń w stępnych. Zaproponowano a lg o ry tm sprzężonych obliczeń pola pręd

kości w ody g ru n to w e j i pola te m p e ra tu ry w zam rażanym gruncie p rzy w y k o rz y s ta n iu m etody elem entów skończonych.

THERMAL ANALYSIS OF U N DERG RO UN D ROCK COOLING PROCESS DU RIN G U N DERG RO UN D WATER MOTION

Sum m ary. The m ethod o f ro ck mass cooling and in it ia l freezing calculations w hen th e v e lo c ity d is trib u tio n o f underg ro u n d w a te r is kn o w n and th e m ethod o f c a lcu la tio n o f th is w a te r ve lo city d is trib u tio n have been described in th e paper. The re s u lts o f these calculations have been presented. A n a lg o rith m o f c a lcu la tio n o f ve lo c ity d is trib u tio n of u n d e rg ro u n d w a te r coupled w it h ca lc u la tio n o f te m p e ra tu re d is trib u tio n d u rin g cooling and in it ia l fre e zin g o f u n d e rg ro u n d rock u sing fin ite elem ent m ethod has been suggested.

THERMISCHE ANALYSE D E S BO D EN SG EFR IER U N G SPR O C ESSES BEI BERÜCKSICHTIGUNG D E R BEW EGUNG DES

UNTERERDISCHEN W ASSERS

Z u sa m m en fa ssu n g. D ie B erechnungsm ethode des Bodensgefrie

rens bei B e rü c k s ic h tig u n g der bekanntes G e sch w in d ig ke itsve rte ilu n g des u n te re rd isch e n W assers un d die B estim m ungsm ethode ihres Ge

schw indigke itsfeldes w u rd e n besprochen. E in ig e Berechnungsergeb

nisse w u rd e n angegeben. D e r A g o rith m u s der B erechnungskopplung des G eschw indig keits- u n d T e m p e ra tu rfe ld e s im gefrierenem Boden bei A nw e n d u n g der M ethode der F in ite n E lem enten w u rd e vorgeschlagen.

(2)

106 Adam Fic, Jan Składzień

1. W STĘP

A n a liz a cieplna procesu ochładzania i zam rażania górotw oru przeprow a

dzana je s t na ogół p rz y założeniu, że w ilgoć z a w a rta w górotworze się nie porusza. T a kie założenie p rz y jm o w a li a u to rzy n in iejszej pracy w dotychczaso

w ych badaniach [1 - 19], podobny przypadek b y ł rów nież ro zp a tryw a n y w pracach in n y c h autorów cytow anych m. in . w [10, 19]. W rzeczywistości woda g ru n to w a znajduje się często w ru ch u . Z w ykle prędkość p rze p ływ u wody je st na ty łe m ała, że nie m a ona istotnego w p ły w u n a proces pow staw ania płaszcza mrożeniowego. W w a rstw a ch g ru n tu o dużej porowatości zdarza się je d n a k, że prędkość ru c h u w ody uzyskuje w artości rzędu 10 m/ds (1 ds = 24 h) lub więcej [21] i wówczas może to w sposób is to tn y negatyw nie w pływ ać na prze

bieg procesu zam rażania górotw oru. W s k ra jn y m p rzyp a d ku może nastąpić sta b iliza cja w y m ia n y ciepła i prze rw a n ie procesu n a ra s ta n ia w a rs tw y lodo- g ru n tu .

W pracy zaproponowano metodę obliczeń zam rażania górotw oru uw zględ

n iającą ko nw ekcyjną w ym ia n ę ciepła zw iązaną z ruchem w ody gruntow ej.

M etodę omówiono w odniesieniu do zadań dw u w ym ia ro w ych dotyczących jednorod nych w a rs tw g ru n tu .

P rzy w y k o n y w a n iu obliczeń liczbow ych rozpatrzono fazę ochładzania oraz wstępnego zam rażania górotw oru. N ie prowadzono obliczeń dla dalszej części procesu, gdy obszar zam rożony m a w y m ia ry znacznie większe n iż średnica o tw o ru mrożeniowego 2 r0. W zw ią zku z ty m było uzasadnione rozpatryw anie prostokątnego obszaru ja k w p. 4. W p rzyp a d ku a n a liz y procesu m rożenia górotw oru związanego z pow staniem klasycznego, zam kniętego płaszcza m ro

żeniowego postępowanie ta k ie nie je s t m ożliw e. Celem niniejszej pracy było je d n a k je d y n ie zbadanie, czy i ja k i w p ły w na zja w iska w y m ia n y ciepła w górotw orze m ają ru c h y wód podziem nych p rz y bardzo m ałej prędkości przesą

czania się cieczy. Do uzyskania tego celu rozpatrzono, ja k ju ż wspomniano, je d y n ie początkow y etap procesu. A n a liz a dalszego przebiegu zam rażania g ó rotw oru w ym aga ju ż zastosowania znacznie bardziej złożonego kodu ob li

czeniowego, uw zględniającego w szystkie r u r y m rożeniowe i obejmującego cały zam rażany i ochładzany obszar. W y n ik a to stąd, że n ie je s t wówczas możliwe w ydzielenie, ja k w [1 - 19], pow tarzalnego, je ś li chodzi o przebieg zjaw isk cieplnych, fra g m e n tu górotw oru.

O bliczenia liczbowe związane z zam rażaniem górotw oru, k tó ry c h re z u lta ty przedstaw iono w niniejszej pracy, prowadzone b y ły w dwóch etapach. W p ie r

w szym określano ro z k ła d prędkości w ody przesączającej się w górotworze.

P rzyjm ow ano p rz y ty m , że pole prędkości cieczy je s t ta k ie ja k w stanie u stalonym , ro zkła d tej prędkości zaś wyznaczano analitycznie . W etapie d ru gim dla znanego pola prędkości w ody określano ro zkła d te m p e ra tu ry w góro

tw orze w yko rzystu ją c metodę elem entów skończonych.

(3)

Analiza cieplna procesu ochładzania górotworu. 107

O ddzielnym fragm entem są num eryczne obliczenia pola prędkości w ody w górotworze d la p rzy p a d k u o p ływ a n ia nieskończenie długiego walca. Również tu w yko rzysta n o metodę elem entów skończonych i p osłużyła ona do w yzna

czenia potencja łu prędkości w ody g ru n to w e j. Przeprowadzone sym ulacje kom puterow e d oprow ad ziły do pow sta n ia koncepcji sprzężonych obliczeń pól pręd

kości w ody g ru n to w e j i te m p e ra tu ry g ó ro tw o ru podczas procesu jego ochładza

nia i m rożenia. P rzedstaw iono zarys te j koncepcji, k tó ra , zgodnie z zam ierze

n ia m i autorów , m a stanow ić podstawę do dalszych badań. B adania ta k ie będą je d n a k w ym a g a ły w y k o rz y s ta n ia k o m p u te ra o znacznie w iększych m ożliw o

ściach obliczeniowych.

2. M E T O D A O B L IC Z E Ń P O L A T E M P E R A T U R Y

R ów nanie przew odzenia ciepła w zam rażanym gruncie zaw ierającym p o ru szającą się wodę g ru n to w ą m a postać:

(pc) + (pc)w v V T = V k VT, (1)

gdzie prędkość v*w o d y g ru n to w e j je s t prędkością pozorną, określoną ja ko stosunek objętości w ody przepływ ającej w jednostce czasu przez pow ierzchnię i pola tej pow ierzchni.

R ów nanie (1) je s t słuszne, p rz y założeniu że te m p e ra tu ra w ody prze p ływ a jącej przez p o ry w gruncie je s t ró w n a te m p e ra tu rz e g ru n tu .

N a brzegach obszaru obliczeniowego będącego w yodrębnio nym fragm en

tem zam rażanego g ó ro tw o ru p rz y jm u je się w a ru n k i brzegowe I I I rodzaju.

W a ru n k ie m początkow ym je s t sta ła te m p e ra tu ra w całym ro zp a tryw a n ym obszarze.

W celu u w zg lę d n ie n ia e fe ktu cieplnego zm ia n y fazy w yko rzystu je się meto

dę uśre d n ie ń całkow ych [14]. Jej is to ta polega n a w p row adzen iu odpowiednio uśrednionej pojem ności cieplnej w łaściw ej w ty c h elem entach, przez któ re przechodzi g ranica m iędzy fazam i. M etoda pozw ala uw zględnić rów nież zróż

nicow anie przewodności cieplnej w obu fazach.

Do przestrzennej d yskre tyza cji zagadnienia zastosowano metodę elemen

tów skończonych z d w u w y m ia ro w y m i iz o p a ra m e tryczn ym i elem entam i czterow ęzłow ym i [14 - 20],

Ze w zględu na n ie w ie lk ie w a rto ści prędkości w w o d y g runtow ej, niezerowy w niezam rożonej części obszaru s k ła d n ik k o n w e kcyjn y w ró w n a n iu (1) tr a k tu je się ta k , ja k w m etodzie elem entów skończonych w ydajność qv wewnętrznego

źródła ciepła. W zw ią zku z ty m w obszarze niezam rożonym p rz y jm u je się

(4)

qv = -(p c)w v V T dla T > T m ,

w obszarze zam rożonym zaś

4v = 0 dla T < T m .

(

2

)

(3) P rzy znanych składow ych v x i v y prędkości v f d la T > T m w dw uw ym iarow ym u kła d zie p ła s k im o trzym u je się

4v = - (pc)v 3T 3 T n

v x T + V , T x 5x y 3y (4)

Pochodne te m p e ra tu ry w tej re la c ji w yznacza się korzystając z in te rp o la cji te m p e ra tu ry stosowanej w m etodzie elem entów skończonych, stąd np.

3T LWy dNj

^ 1 d x ’ d x i = l

(5)

3. M E T O D A O B L IC Z E Ń P O L A P R Ę D K O Ś C I W O D Y G R U N T O W E J Rozkład prędkości pozornej U*ruchu w ody g ru n to w e j w ośrodku zia rn istym , ja k im je s t g ru n t (górotw ór), opisany je s t ró w n a n ie m ciągłości s tru g i i rów na

niem zachowania pędu. W p rzyp a d ku ustalonego prze p ływ u p ły n u n ieściśli

wego ró w n a n ie ciągłości s tru g i m a postać [22]

V ^ = 0 , (6)

ró w n a n ie zachowania pędu n a to m ia s t prow adzi do zależności wyrażającej praw o D arcy’ego

v*= - k Vh, (7)

gdzie wysokość rozporządzalna w ody w głębnej je s t określona wzorem

h = z + — . (8)

P wg

S kojarzenie re la c ji (6) i (7) daje ró w n a n ie opisujące ro zkła d wysokości rozporządzalnej h

Vk Vh = 0 . (9)

(5)

Wzór (7) stanow i podstawę do fo rm u ło w a n ia w a ru n k ó w brzegowych dla ró w nania (9) oraz do w yznaczania pola prędkości wody, gdy znany je s t ro zkła d wysokości rozporządzalnej h.

Jeżeli ośrodek je s t jednorod ny, to ró w n a n ie (9) m ożna zapisać w form ie

W a ru n k i brzegowe zw ykle fo rm u łu je się w yko rzystu ją c znajomość składo

wej n o rm a ln e j prędkości w ody na brzegu rozpatryw ane go obszaru. N a podsta

wie (12) składow a ta w y n ik a z ró w n a n ia

którego postać je s t analogiczna, ja k ró w n a n ia w a ru n k u brzegowego I I rodzaju dla ró w n a n ia przew odzenia ciepła.

Zagadnienie brzegowe (10), (13) nie m a jednoznacznego rozw iązania za względu na potencjał (p. Do w yznaczenia ro z k ła d u prędkości w ody w ystarcza je d n a k ty lk o znajomość ro z k ła d u p a ra m e tru tp z dokładnością do stałej. W ob

liczeniach num erycznych najdogodniej je s t p rzyją ć pew ną w artość fu n k c ji w w yb ra n ym p unkcie s ia tk i.

Zagadnienie brzegowe (10), (13) opisane je s t ró w n a n ia m i analogicznym i do rów nań opisujących pole te m p e ra tu ry p rz y w ystę p o w a n iu typow ych w a ru n ków brzegowych. Może być więc ono rozw iązyw ane za pomocą m etody elemen

tów skończonych, podobnie ja k wyznaczane są pola te m p e ra tu ry i p rzy zasto

sow aniu np. p ro g ra m u w ykorzystyw anego do obliczeń dotyczących procesu zam rażania g ó rotw oru m etodą opisaną w p. 2.

Znając num eryczne rozw iązanie w o dniesie niu do fu n k c ji cp można łatw o wyznaczyć składowe vx i v y prędkości w ody g ru n to w e j w dowolnym punkcie.

W ty m celu w yko rzystu je się praw o D a rcy’cgo oraz stosowaną w metodzie elem entów skończonych in te rp o la cję fu n k c ji <p, skąd w y n ik a np.:

V 2cp=0,

(

₁₀

)

gdzie (p je s t potencjałem prędkości wód w głębnych

(p = k h ,

(

11

)

wzór (7) zaś sprowadza się do postaci

(

12

)

n ■ V(p + v n = 0 , (13)

(14)

(6)

4. R E Z U LT A T Y W S T Ę P N Y C H O B L IC Z E Ń P O L A T E M P E R A T U R Y I P O LA P R Ę D K O Ś C I W O D Y W GÓROTW ORZE

W c e lu z ilu stro w a n ia problem u w ykonano przykładow e obliczenia dla pro

cesu zam rażania g ó rotw oru w okół pojedynczej r u r y m rożeniowej w jednorod

nej w a rs tw ie , w k tó re j przepływ a woda g ru n to w a z prędkością v „ w dostatecz

n ie dużej odległości od r u r y m rożeniowej i zamrożonej w a rs tw y [20]. P rzyję ty obszar obliczeniow y pokazano na rys. 1. Obszar te n sta n o w i połowę zaizolowa

nego cieplnie k w a d ra tu o rozm iarach 2 x 2 m. W jego środku znajduje się ru ra m rożeniow a o p ro m ie n iu zew nętrznym r 0 = 0,13 m, przez k tó rą w c h w ili początkowej ( t = 0) zaczyna przepływ ać chłodziw o o te m p e ra tu rze -30°C , p rzy czym w spółczynnik p rz e n ik a n ia ciepła od g ó rotw oru do cieczy chłodzącej w ynosi 500 W /(m ■ K). P rzyjęto te m p e ra tu rę początkową górotw oru 20°C, objętościową e n ta lp ię w ła ściw ą zm ia n y fazy 1,3 • 108 J /m 3, w spółczynnik prze

w odzenia ciepła d la góro tw o ru zamrożonego i niezamrożonego odpowiednio 2,2 i 1,2 W /(m • K ), a objętościową pojemność cieplną w łaściw ą w tych w a ru n kach 3,2 106 i 2,2 106 J/(m 3 K). O bliczenia w ykonano dla dwóch w a ria n tów , p rzyjm u ją c b ra k ru c h u w ody podziem nej oraz = 40 m/ds. Rozkład prędkości w ody g runtow ej w o kó ł r u r y lu b zamrożonej w a rs tw y przyjęto ja k p rz y bezcyrku la cyjn ym opływ ie w ośrodku nieskończonym [23] w alca kołowe

go o p ro m ie n iu ró w n ym uśrednionem u p ro m ie n io w i zamrożonej w a rstw y.

N a rys. 2 pokazano otrzym ane p rzykładow o granice zamrożonej w a rs tw y w obu w a ria n ta c h d la w yb ra n ych czasów. W w a ria n cie z poruszającą się wodą

Rys. 1. Podział rozpatrywanego fragmentu górotworu na elementy skończone Fig. 1. Subdivision of the examined fragment of underground rock into finite elements

(7)

gruntow ą obszar zam rożony m a k s z ta łt zbliżo n y do e lip sy o środku przesunię

tym w stosunku do środka r u r y m rożeniow ej w k ie ru n k u zgodnym z k ie ru n kiem dopływającej wody. R e z u lta t te n je s t zgodny z o czekiw ania m i i z in fo r

m acjam i podan ym i w [21].

■i - i = 144 h

Rys. 2. Przykładowa postać granicy obszaru zamrożonego (1 ds = 24 h; współczynnik przejmowania ciepła)

Fig. 2. Exemplary form of the frozen area boundary (1 ds = 24 h ; - heat transfer coefficient)

Za pomocą tego samego program u, w ykorzystującego metodę elementów skończonych do w yznaczania p o te n cja łu prędkości w ody g runtow ej cp i p rzy

stosowanego do obliczeń na tej podstaw ie pola prędkości, w ykonano w analo

gicznym obszarze obliczenia ro z k ła d u prędkości w ody g runtow ej. Przyjęto, że na brzegu A B , CD i E F składow a n o rm a ln a prędkości v n je s t ró w n a zeru, a na brzegu A F w ynosi -v„>. W a ru n k i brzegowe mogą być określone w ta k i sposób w p rzyp a d ku o p ływ a n ia przez wodę g ru n to w ą rzędu r u r m rożeniow ych lub otaczających te r u r y w alcow ych w a rs tw lodogrun tow ych, je ż e li woda g ru n to wa dopływ a prostopadle do tego rzędu. Obszar obliczeniow y je s t p o w ta rza l

nym fragm entem niezam rożonej części g ó ro tw o ru otaczającego ru rę . Oblicze

nia w ykonano dla różnych p ro m ie n i r 0 półkolisteg o brzegu BC. N ie któ re re z u lta ty tych obliczeń w ra z z prędkościam i v 0 w yznaczonym i dla bezcyrku- lacyjnego opływ u w alca w ośrodku nieskończonym [23] zaw arto w tablicach 1 i 2. P orów nanie otrzym a n ych re z u lta tó w p o tw ie rd za poprawność zastosowa

nej m etody (dobra zgodność w y n ik ó w w p o b liżu brzegu BC). Znaczna ich rozbieżność w pewnej odległości od g ra n ic y BC, szczególnie dla w iększych

(8)

je s t w y k o rzystyw a n ie wspom nianego analitycznego rozw iązania opisującego ro z k ła d prędkości, je ż e li sum a p ro m ie n ia r u r y i grubości płaszcza mrożenio- wego staje się p o rów nyw aln a z odległością pom iędzy ru ra m i m rożeniow ym i.

T a b lic a 1 S k ła d o w e p r ę d k o ś c i w o d y g r u n to w e j w ś r o d k a c h e le m e n tó w d la ro = 0,13 m,

w m /d s (1 d s = 24 h)

element vx vy ^{V x 0} ^{V y 0}

1 39,85 0,4041 39,31 0,3673

2 39,45 0,5355 38,88 0,5871

3 38,80 0,8478 38,18 0,9441

4 37,69 1,394 37,00 1,530

5 35,77 2,323 34,99 2,502

6 32,35 3,926 31,44 4,141

7 26,07 6,704 24,98 6,956

8 14,07 11,31 12,60 11,90

9 40,16 0,4406 39,85 0,5138

10 40,19 0,7209 39,75 0,8392

11 40,07 1,261 39,58 1,388

12 39,79 2,174 39,28 2,333

13 39,26 3,773 38,72 4,010

14 38,26 6,709 37,62 7,106

15 36,24 12,44 35,33 13,14

16 31,36 24,45 39,94 25,89

17 40,51 0,2203 40,15 0,5138

18 40,75 0,6611 40,25 0,8392

19 40,98 1,247 40,42 1,388

20 41,31 2,180 40,72 2,333

21 41,87 3,793 41,28 4,010

22 42,94 6,751 42,38 7,106

23 45,09 12,54 44,67 13,14

24 50,34 25,32 50,06 25,89

25 41,40 0,1943 40,69 0,3673

26 41,70 0,4664 41,12 0,5871

27 42,32 0,8288 41,82 0,9441

28 43,41 1,394 43,00 1,530

29 45,31 2,335 45,01 2,502

30 48,69 3,957 48,56 4,141

31 54,87 6,831 55,02 6,956

32 66,55 11,68 67,40 11,90

33 41,40 -0,1943 40,69 -0,3673

34 41,70 -0,4664 41,12 -0,5871

35 42,32 -0,8288 41,82 -0,9941

36 43,41 -1,394 43,00 -1,530

37 45,31 -2,335 45,01 -2,502

38 48,69 -3,957 48,56 -i,141

39 54,87 -6,831 55,02 -6,956

40 66,56 -11,68 67,40 -11,90

(9)

c d . t a b l i c y 1

element vx vy Vx0 vy0

41 40,51 -0,2203 40,15 -0,5138

42 40,75 -0,6611 40,25 -0,8392

43 40,98 -1,247 40,42 -1,388

44 41,31 -2,180 40,72 -2,333

45 41,87 -3,793 41,28 -4,010

46 42,94 -6,751 42,38 -7,106

47 45,09 -12,54 44,67 -13,14

48 50,34 -25,32 50,06 -25,89

49 40,16 -0,4406 39,85 -0,5138

50 40,19 -0,7209 39,75 -0,8392

51 40,07 -1,261 39,58 -1,388

52 39,79 -2,174 39,28 -2,333

53 39,26 -3,773 38,72 -4,010

54 38,26 -6,709 37,62 -7,106

55 36,24 -12,44 35,33 -13,14

56 31,36 -24,45 29,94 -25,89

57 39,85 -0,4041 39,31 -0,3673

58 39,45 -0,5355 38,88 -0,5871

59 38,80 -0,8478 38,18 -0,9441

60 37,69 -1,394 37,00 -1,530

61 35,77 -2,323 34,99 -2,502

62 32,35 -3,926 31,44 -4,141

63 26,07 -6,704 24,98 -6,956

64 14,07 -11,31 12,60 -11,90

T a b l i c a 2 S k ła d o w e p r ę d k o ś c i w o d y g r u n t o w e j w ś r o d k a c h e l e m e n t ó w d l a ro = 0,5 m ,

w m /d s (1 d s = 24 h )

element vx vy Vx0 vy0

1 38,87 7,492 30,99 4,763

2 36,85 7,697 29,33 5,533

3 34,40 8,322 27,30 6,447

4 31,40 9,347 24,77 7,538

5 27,65 10,76 21,61 8,844

6 22,89 12,57 17,61 10,41

7 16,79 14,74 12,46 12,31

8 8,919 17,17 5,762 14,59

9 42,44 6,080 38,01 6,733

10 43,47 7,752 37,53 8,152

11 43,92 9,932 36,88 9,961

12 43,84 12,77 36,01 12,30

13 43,18 16,48 34,81 15,37

14 41,85 21,42 33,11 19,49

15 39,57 28,11 30,61 25,14

16 35,68 37,42 26,80 33,10

17 48,05 1,992 41,99 6,733

(10)

114 Adam Fie, Jan Sktadzien

c d . ta b l i c y 2

element v* vy VxO VyO

18 50,27 4,903 42,47 8,152

19 52,18 8,037 43,12 9,961

20 53,97 11,61 43,99 12,30

21 55,81 15,91 45,19 15,37

22 57,95 21,37 46,89 19,49

23 60,75 28,71 49,39 25,14

24 64,80 39,26 53,20 33,10

25 65,34 2,420 49,01 4,763

26 66,36 4,162 50,67 5,533

27 68,00 5,962 52,70 6,447

28 70,36 7,895 55,23 7,538

29 73,59 10,04 58,39 8,844

30 77,90 12,45 62,39 10,41

31 83,59 15,17 67,54 12,31

32 91,05 18,06 74,24 14,59

33 65,34 -2,420 49,01 -4,763

34 66,36 -4,161 50,67 -5,533

35 68,00 -5,962 52,70 -6,447

36 70,36 -7,895 55,23 -7,538

37 73,59 -10,03 58,39 -8,844

38 77,90 -12,45 62,39 -10,41

39 83,59 -15,17 67,54 -12,31

40 91,05 -18,06 74,24 -14,59

41 48,05 -1,992 41,99 -6,733

42 50,27 -4,903 42,47 -8,152

43 52,18 -8,037 43,12 -9,961

44 53,97 -11,61 43,99 -12,30

45 55,81 -15,91 45,19 -15,37

46 57,95 -21,37 46,89 -19,49

47 60,75 -28,71 49,39 -25,14

48 64,80 -39,26 53,20 -33,10

49 42,44 -6,080 38,01 -6,733

50 43,47 -7,752 37,53 -8,152

51 43,92 -9,932 36,88 -9,961

52 43,84 -12,77 36,01 -12,30

53 43,18 -16,48 34,81 -15,37

54 41,85 -21,42 33,11 -19,49

55 39,57 -28,11 30,61 -25,14

56 35,69 -37,42 26,80 -33,10

57 38,87 -7,492 30,99 -4,763

58 36,84 -7,697 29,33 -5,533

59 34,40 -8,322 27,30 -6,447

60 31,40 -9,347 24,77 -7,538

61 27,65 -10,76 21,61 -8,844

62 22,89 -12,57 17,61 -10,41

63 16,79 -14,73 12,46 -12,31

64 8,919 -17,17 5,762 -14,59

(11)

Rys. 3. Przykładowa postać elementów siatki wtórnej stanowiących fragmenty elementów siatki pierwotnej

Fig. 3. Exemplary form of second net elements which are fragments of primary net elements

5. SPRZĘŻONE O B L IC Z E N IA P O L A P R Ę D K O Ś C I W O D Y I PO LA T E M P E R A T U R Y W GÓROTW ORZE

P rzy zastosow aniu zaproponow anych m etod można, ja k w ykazano, w yko rzystywać typow e p rocedury m etody elem entów skończonych do rozw iązyw a

nia zadań przew odzenia ciepła w celu ro zw ią za n ia zagadnienia m rożenia górotw oru w p rz y p a d k u w ystę p o w a n ia ru c h u wód podziem nych. Pewne pro

blem y p o ja w ia ją się je d n a k podczas p ra ktyczn e j re a liz a c ji a lg o ry tm u obliczeń.

M ożna założyć, że szybkość n a ra s ta n ia płaszcza m rożeni owego je s t na ty le mała, iż ro z k ła d prędkości w ody g ru n to w e j w każdej c h w ili m ożna przyjąć jako ro z k ła d u s ta lo n y odpow iadający opływ ow i a k tu a ln e j w a rs tw y lodogrun- towej. Wówczas w ko le jn ych k ro k a c h czasu m ożna wyznaczać n a jp ie rw roz

k ła d prędkości w ody g ru n to w e j om ywającej ru rę m rożeniow ą lu b zamrożoną w arstw ę, a następn ie ro z k ła d te m p e ra tu ry w górotw orze. R ozkład te m p e ra tu ry w yznacza się je d n a k zarówno w zam rożonej, ja k i niezam rożonej części g ru n tu , n a to m ia s t ro z k ła d prędkości je d y n ie w części niezam rożonej.

W elem entach, w k tó ry c h przebiega granica zm ia n y fazy, m ożna zaniedbać składow ą unoszenia w ró w n a n iu przew odzenia ciepła. Prędkość w ody g ru n to wej je s t bow iem w otoczeniu g ra n icy zm ia n y fazy styczna do tej granicy, tzn.

do iz o te rm y krzepnięcia, skąd v V T = 0 w w w . elem entach.

W dążeniu do m in im a liz a c ji czasu obliczeń ko rzystn e je s t w yznaczanie pola prędkości w ody g ru n to w e j p rz y założeniu postaci s ia tk i w tó rn e j m ożliw ie n ajbardzie j zbliżonej do postaci s ia tk i p ie rw o tn e j, p rzyję te j p rz y obliczaniu rozkładów te m p e ra tu ry w k o le jn y c h k ro k a c h czasu. Proponuje się, aby sia tka w tó rn a sk ła d a ła się z niezam rożonych elem entów s ia tk i p ierw otnej oraz z elem entów reprezentu jących niezam rożone części ty c h elem entów s ia tk i p ie r

w otnej, przez k tó re przebiega g ranica m iędzy fazam i. P rzy stosow aniu ele

(12)

m entów izoparam etrycznych czterowęzłowych m o żliw a je s t reprezentacja w w . części ta k ic h elem entów przez elem ent tró jk ą tn y , czw orokątny lu b ele

m e n t tró jk ą tn y i czw orokątny (elem enty e’ i e” na rys. 3).

Każdorazowe obliczenia ro z k ła d u prędkości p o w in n y być poprzedzone zde

fin io w a n ie m s ia tk i w tó rn e j, co polega na o kre śle n iu n u m e ra cji elementów, w ęzłów i połączeń w ęzłów w tej siatce. Konieczne je s t też przyporządkow anie w ęzłom s ia tk i p ie rw o tn e j w ęzłów s ia tk i w tó rn e j. Stosowna procedura może polegać na a n a lizo w a n iu kole jn ych elem entów s ia tk i pierw otnej i w artości te m p e ra tu ry w węzłach. Jeżeli przez elem ent przebiega granica zm iany fazy, to m a miejsce określenie je j przebiegu oraz w yznaczenie dodatkowego, w p orów n a n iu z sia tk ą pie rw o tn ą , elem entu (lu b elem entów - rys. 3) i dodatko

w ych węzłów. E le m e n ty w tórne, utw orzone z niezam rożonych części elemen

tó w p ie rw o tn ych , przez k tó re przebiega g ra n ica zm ia n y fazy, m ają w pływ je d y n ie na składow ą objętościową m acierzy w spółczynników . E le m e n ty s ia tk i p ie rw o tn e j p o w in n y być num erow ane ta k , aby kolejno generowane węzły s ia tk i w tó rn e j g w a ra n to w a ły m in im a ln ą szerokość pasm a m acierzy współ

czynników . W a ru n e k te n spełnia n um eracja pokazana na rys. 1.

Rozkład prędkości w ody g ru n to w e j m ożna wyznaczać i korygować po u p ły w ie pewnej liczb y k ro k ó w czasowych. W te n sposób je s t m ożliw e skrócenie czasu obliczeń p rz y p ra ktyczn ie niezm iennej ich dokładności.

O Z N A C Z E N IA

g - przyśpieszenie ziem skie,

h - wysokość rozporządzalna w ody w głębnej, k - w sp ó łczyn n ik przepuszczalności,

L W — liczba w ęzłów w elemencie,

ń* - w ersor norm a ln e j zew nętrznej do brzegu, N j - fu n k c ja k s z ta łtu zw iązana z i-ty m węzłem,

p - ciśnienie,

< ł v - wydajność w e w n ę trzn ych źródeł ciepła,

T, Tm, Ti - te m p e ra tu ra , te m p e ra tu ra z m ia n y fa z y i te m p e ra tu ra w i-ty m węźle,

vx> vy _ prędkość pozorna w ody i je j składowe w dwuw ym iarow ym układzie p ła skim ,

x, y - w spółrzędne w u kła d zie p ła skim , z - w spółrzędna pionowa,

V - operator N abla,

(p - potencjał prędkości wód w głębnych, X — w sp ó łczyn n ik przew odzenia ciepła, P> Pw - gęstość i gęstość w ody g runtow ej,

(pc) — objętościowa pojemność cieplna w łaściw a,

(13)

(pc)„, - objętościowa pojemność cieplna w ła ściw a wody, x - czas.

L IT E R A T U R A

[1] S kładzień J.: Zastosowanie m etody różnicow ej do a n a lizy cieplnej proce

su zam rażania górotw oru, R e fe ra ty S ym pozjum W y m ia n y C iepła i M a sy, W arszaw a - Jabłonn a 1976.

[2] C ie rp ka W ., S kładzień J.: Prognozowanie procesu zam rażania górotw oru do dużych głębokości za pomocą m etod num erycznych - różnicowych, R eferaty K on fe re n cji „G łębokie m rożenie g ó ro tw o ru p rz y głębieniu szy

bów” , Częstochowa 1976.

[3] S kładzień J.: Zastosowanie m etody różnicow ej do a n a liz y cieplnej proce

su zam rażania g ó ro tw o ru za pomocą podwójnego k rę g u otworów, Zeszy

ty N aukow e Pol. Śl., E n e rg e tyka z. 63, G liw ice 1978.

[4] S kładzień J.: A n a liz a procesu za m ra ża n ia g ó ro tw o ru podw ójnym k rę giem otw orów m rożeniow ych, R e fe ra ty X Zjazdu T erm odynam ików , W rocław 1978.

[5] S kładzień J.: A n a liz a cieplna procesu zam ra ża n ia górotw oru p rz y pomo

cy m etody różnicow ej, A rc h iw u m T e rm o d y n a m ik i i Spalania 1, 1978.

[6] S kładzień J.: B ezw ym iarow a a n a liza cieplna procesu zam rażania góro

tw o ru , Zeszyty N aukow e Pol. Śl., E n e rg e tyka z. 67, G liw ice 1978.

[7] S kładzień J.: Z am rażanie g ó ro tw o ru p rz y ograniczonej mocy agregatów m rożeniow ych, Zeszyty N aukow e Pol. Śl., E n e rg e tyka z. 71, G liw ice 1979.

[8] S kładzień J.: The a p p lica tio n o f th e difference equations in th e rm a l analysis o f rock-freezing, re fe ra t w ygłoszony n a k o n fe re n cji „N u m e rica l M ethods in T h e rm a l P roblem s” , Swansea 1979.

[9] S kładzień J . : M etodyka o kre śla n ia o p tym a ln ych p a ra m e tró w m rożenia górotw oru, M a te ria ły Zjazdowe X I Zjazdu T erm odynam ików , Szczecin - Ś w inoujście 1981.

[10] S kładzień J.: A n a liz a cieplna i ekonom iczna m rożenia górotw oru, Zeszy

t y N aukow e Pol. Śl., E n e rg e tyka z. 78 (m onografia), G liw ice 1981.

[11] S kładzień J.: Dobór mocy zespołu agregatów do głębokiego zam rażania górotw oru, M a te ria ły K o n fe re n cji N aukow o -T e ch n iczn e j „K o m p u te ry w P ro je kto w a n iu i E ksp lo a ta cji U rządzeń C hłodniczych i K lim a ty z a c y j

n ych” , X V I I I D n i C hłodnictw a , Poznań 1984.

[12] S kładzień J.: A n a liz a cieplna m rożenia górotw oru, B udow nictw o W ęglo

we - P ro je k ty - P roblem y 5, 1984.

[13] S kładzień J.: B ezw ym iarow a a n a liza cieplna m rożenia górotw oru je d n ym kręgiem ru r, B u dow nictw o W ęglowe - P ro je k ty - Problem y 2,1987.

(14)

[14] F ic A.: B łę d y m odelow ania cieplnego e fe ktu zm ia n y fazy p rz y stosowa

n iu m etody uśre d n ie ń całkow ych do ro zw ią zyw a n ia zagadnień Stefana, Z biór R eferatów X X V I Sym pozjonu „M odelow anie w M echanice” , G liw i

ce - W isła 1987.

[15] Fic A., S kładzień J.: O bliczanie zam rażania górotw oru p rzy w yko rzysta n iu m etody u średnień całkowych, Z b ió r R eferatów X X V III Sympozjonu

“M odelow anie w M echanice” , G liw ice — W is ła 1989.

[16] F ic A., S kładzień J.: A n a liz a num eryczna zam rażania górotw oru przy stosow aniu różnych u kła d ó w r u r m rożeniow ych, R eferaty V I I Sympo

zjum W y m ia n y C iepła i M asy, W arszaw a — J a d w is in 1989.

[17] F ic A.: W ybrane aspekty stosowania m etody u średnień całkow ych do ro zw ią zyw a n ia zagadnień Stefana, Zeszyty Naukow e Pol. SI., E nergety

k a z. 107, G liw ice 1989.

[18] Fic A., S kładzień J.: W p ły w u k ła d u r u r na przebieg zam rażania górotwo

r u podw ójnym kręgiem , M a te ria ły K onferencyjne X IV Zjazdu Term ody

na m ikó w , K ra k ó w 1990.

[19] S kładzień J., F ic A.: Zastosowanie M E S do a n a lizy cieplnej procesu zam rażania g ó rotw oru p rz y różnej k o n fig u ra c ji r u r m rożeniowych, Ze

szyty N aukow e Pol. Śl., E nergetyka, G liw ice (przyjęto do d ru k u w 1994 r.).

[20] Fic A., S kładzień J.: Koncepcja obliczeń zam rażania górotw oru z uw z

ględnieniem ru c h u wód podziem nych, R eferaty V I I I Sym pozjum W ym ia n y C iepła i M asy, B iałow ieża 1992.

[21] T ru p a k N. G.: Z a m orażiw anije g ru n tó w p r i s tro itie ls tw ie podziemnych soorużenij, N ie d ra , M oskw a 1979.

[22] Eagleson P. S.: H yd ro lo g ia dynam iczna, PW N, W arszaw a 1978.

[23] Gryboś R.: Podstaw y m e ch a n iki płynów , PW N , W arszaw a 1989.

Recenzent: D r inż. A n to n i G uzik

W p łynęło do R edakcji: 31. 10. 1994 r.

A b stract

The underg ro u n d rock cooling and freezing process is u s u a lly a subject of th e rm a l analysis i f th e re are n o t u n d e rg ro u n d w a te r movements. Such as

su m p tio n was accepted by th e authors in th e ir previous papers [1 - 19]. In the papers [1 - 13] th e fin ite difference equations o f elem e n ta ry balance w ith the fo rw a rd difference a p p ro xim a tio n were used and th e underground rock freez

in g process w ith single and double circle o f freezing pipes was considered. In

(15)

the papers [14 - 19] th e fin ite elem ent m ethod w ith th e in te g ra l averaging method fo r ta k in g in to account th e la te n t h e a t was applied and the u n d e r

ground ro ck freezing process w ith single and double circle o f freezing pipes was also considered. In th e paper [20] th e suggestion o f th e rm a l calculations for undergro und rock cooling and freezing process w ith w a te r flo w was presented. T h is suggestion is here a subject o f considerations togethe r w ith some re su lts o f exem plary, in tro d u c to ry calculation s. A lso th e fin ite elem ent method w ith th e in te g ra l a veraging m ethod fo r ta k in g in to account th e la te n t heat was applied. F ig 1 shows th e s u b d ivisio n o f h a lf o f the area beeing analysed in to fin ite elem ents. O n ly th e f ir s t stage o f th e undergro und rock cooling and freezing process is considered. T h is is w h y th e analysed space in Fig. 1 has a h a lf-s q u a re fo rm w ith a h a lf o f fre e zin g pipe in th e m id d le o f the inner longer side. Fig. 2 presents th e e xe m p la ry fo rm o f th e frozen area boundary fo r 2 cases: i f th e re are n o t undergro und w a te r movem ents and when the in it ia l v e lo c ity o f th e w a te r is v« = 40 m/ds (1 ds = 24 h). The calculation o f w a te r ve lo city d is trib u tio n is th e second problem w h ich is considered in th e paper. I t was proved th a t th e set o f equations describing w ater m ovem ents has th e analog ical fo rm to th e set o f equations fo rm in g the therm al b o u n d a ry problem . T h is fa c t enables one to use th e same procedures as before w ith th e fin ite elem ent m ethod. The re s u lts o f exem plary calcula

tions are given in Tab. 1 and 2. These re s u lts have a fo rm o f values o f ve locity components in every elem ent. In Tab. 1 and 2 co m p a ra tive ly also th e results of approxim ate a n a ly tic calculation s fo r cyrcu la tio n le ss w a te r flo w ro u n d the cylinder in in fin ite space are presented. In th e la s t chapter an a lg o rith m of calculation o f ve lo city d is trib u tio n o f u n d e rg ro u n d w a te r coupled w ith calcula

tion o f te m p e ra tu re d is trib u tio n d u rin g cooling and in it ia l freezing o f u n d e r

ground rock u sin g fin ite elem ent m ethod has been suggested.