MODELOWANIE  MATEMATYCZNE W GEOLOGII

Andrzej Pawuła

Autor  odbył staż naukowy w zakresie modelowania matematycznego w hydrogeologii w Ecole Nationale Superiere des Mines de Paris, Centrum Informatyki Geologicznej w Fontainebleau  (1976) i wykładał hydrodynamikę wód podziemnych oraz modelowanie matematyczne na Uniwersytecie w Constantine w latach 1981-1991, jest   autorem kilku  programów komputerowych, m.in. modelu impulsacyjnego zlewni hydrologicznej, którego wynikiem jest hydrogram odpływu oraz programu  obliczeniowego stateczności skarpy.  W ramach prac badawczych Instytutu Kształtowania Środowiska (1972-1977) autor kierował zespołem, który  opracował regionalny model symulacyjny systemu geohydraulicznego Środkowej Wielkopolski. 


MODEL MATEMATYCZNY POLA FILTRACJI

Treść: Hipoteza Darcy, współczynnik filtracji, twierdzenie Bernouliego,  ogólne równanie filtracji w nieustalonym reżimie przepływu - dla warstwy pod ciśnieniem i o zwieciadle swobodnym, rozwiązanie numeryczne równania filtracji dla przepływu ustalonego w układzie dwuwymiarowym (metoda różnic skończonych), siatka dyskretyzacyjna pola filtracji i rozwiązanie układu równań (metoda iteracyjna).      

Założenie wyjściowe:  hipoteza Darcy o ciągłości ośrodka wodonośnego i fikcyjna prędkość filtracji ():

        (5.1)

gdzie: k - współczynnik filtracji (m/s); H - wysokość hydrauliczna (potencjał hydrauliczny) (m)

Równanie ciągłości przepływu masy, w którym gęstość masowa wody (ro) jest wielkością skalarną:

                (5.2)

Zgodnie z twierdzeniem Bernouliego wysokość hydrauliczna (H) jest funkcją wysokości położenia (z), wysokości ciśnienia (p/ ro g) i wysokości prędkości (v /2g). W ośrodku gruntowym przy względnie małych prędkościach przepływu wysokość prędkości jest zaniedbywalnie mała, tak więc wysokość hydrauliczna zdefiniowana jest formułą:

         (5.3)

gdzie: z - wysokość położenia (m); p - ciśnienie (N m-2 ); ro  - gęstość masowa wody (kg/m3 lub N s2 m-4 ); g - przyspieszenie siły ciężkości (m s-2 ); u - prędkość przepływu wody (m s-1)

Równanie różniczkowe cząstkowe (5.2) oraz prawo Darcy (5.1) stanowią wyjściowy opis matematyczny filtracji wody jako przepływu jednofazowego, w reżimie nieustalonym. W celu zastosowania powyższego układu równań do obliczania przepływu w warstwie wodonośnej o określonej miąższości (m), korzysta się z relacji charkterystycznej dla ośrodka elastycznego oraz z definicji współczynnika zasobności sprężystej warstwy wodonośnej. Ośrodek elastyczny charakteryzują dwie definicje współczynnika ściśliwości (beta), który jest odwrotnością masowego modulu elastyczności:

             (5.4)

gdzie: ro  - gęstość ośrodka; V - objętość ośrodka; p - ciśnienie

Z porównania powyższych równań  (5.4) wynika, że zmiany objętości ośrodka (V) są proporcjonalne do zmian jego gęstości (ro), lecz ze znakiem przeciwnym:

                       (5.5)

Współczynnik zasobności sprężystej warstwy wodonośnej (S), definiuje się natomiast jako stosunek wody uwolnionej lub zmagazynowanej, przez jednostkę powierzchni wodonośca (P) wyrażoną w m2 , do odpowiadającej zmiany potencjału hydraulicznego:

                                 (5.6)

Ponieważ iloczyn powierzchni warstwy wodonośnej (P) i jej miąższości (m) stanowi objętość warstwy wodonośnej P m = V, stąd P = V/m, można przekształcić formułę (5.6) do postaci:

                                 (5.7)

i uzyskane wyrażenie    wprowadzić do równania (5.5), które przyjmuje postać:

                             (5.8)

Stąd otrzymuje się wyrażenie , które można wprowadzić do równania (5.2), otrzymując:

         (5.9)

Przyjmując założenie upraszczające, że woda jest nieściśliwa czyli zakładając, że gęstość masowa wody jest wielkością stałą (ro = constans), można ten parametr wyciągnąć przed znak dywergencji i równanie filtracji uprościć do postaci:

                         (5.10)

Podstawiając formułę Darcy (5.1) i zakładając izotropowość współczynnika filtracji oraz upraszczając znaki, otrzymuje się postać ogólną równania filtracji dla przepływu nieustalonego:

         (5.11)

Parametrami w tym równaniu jest współczynnik filtracji (k), miąższość warstwy wodonośnej (m) oraz współczynnik zasobności sprężystej (S). Wielkością poszukiwaną jest wysokość hydrauliczna (H). Ogólne równanie filtracji (5.11) wymaga adaptacji dla określonych warunków hydrogeologicznych.

Poziom swobodny:

W przypadku wodonośca o zwierciadle swobodnym miąższość warstwy (m) należy traktować jako wielkość zmienną od wysokości hydraulicznej (H), dlatego parametr ten wymaga wstawienia pod znak dywergencji:

    (5.12)

Poziom ciśnieniowy:

W przypadku poziomu ciśnieniowego można założyć, że miąższość warstwy wodonośnej (m) jest wielkością stałą, podobnie jak współczynnik filtracji (k) i można przedstawić ich iloczyn w postaci współczynnika wodoprzewodności (T = k m), wyłączonego przed znak dywergencji i otrzymujemy równania:

            (5.13)

Reżim nieustalony:

W układzie dwuwymiarowym równanie filtracji dla przepływu nieustalonego w warstwie wodonośnej, redukuje się do wzoru:

        (5.14)

gdzie: Tx,Ty - współczynniki wodoprzewodności, uwzględniające heterotropowy charakter ośrodka wodonośnego.

Reżim ustalony:

Równanie na przepływ w reżimie ustalonym otrzymuje się przez uproszczenie równania (5.12):

                          (5.15)

Wprowadzając analogiczne przekształcenia otrzymuje się dla poziomej warstwy wodonośnej o miąższości (m) równanie filtracji typu parabolicznego (równanie Laplace'a):

                        (5.16)

Dla układu dwuwymiarowego, w przekroju pionowym - w miejsce miąższości warstwy wodonośnej występuje jednostkowa szerokość strumienia filtracji l = 1 m :

                                (5.17)

System wielowarstwowy:

Problem modelowania systemów wielowarstwowych rozwiązywany jest przez rozwinięcie wyrazu q w powyższym równaniu, a który określa wydajność przepływu między elementarną objętością rozpatrywanego układu a układem zewnetrznym. Może to być przepływ założony, oznaczający pobór wody z warstwy wodonośnej lub iniekcję wody, a może to również być przepływ obliczany, uwarunkowany różnicą ciśnień i oporem filtracyjnym warstwy rozdzielającej. W tym przypadku wydajność przepływu q można określić wzorem:

                                (5.18)

gdzie: qe - wydajność poboru lub iniekcji wody do warstwy wodonośnej; qi - wydajność infiltracji; qp - wydajność przesączania między przedmiotową warstwą a sąsiednią warstwą wodonośną.

Przesączanie wody przez warstwę półprzepuszczalną można obliczyć wzorem:

                                (5.19)

gdzie: Tp - współczynnik wodoprzewodności warstwy półprzepuszczalnej; H - wysokość hydrauliczna w punkcie x,y rozpatrywanej warstwy; Ho - wysokość hydrauliczna w sąsiedniej warstwie wodonośnej

Przy tak sformułowanym równaniu (5.19) przepływ zasilający przedmiotową warstwę wodonośną będzie miał wartość dodatnią. Współczynnik wodoprzewodności warstwy półprzepuszczalnej obliczany jest wzorem:

        (5.20)

gdzie: ks - współczynnik filtracji warstwy półprzepuszczalnej (m/s); ms - miąższość warstwy półprzepuszczalnej (m); P - powierzchnia przekroju jednostkowego (m )

c.d. w opracowaniu


Wyjaśnienia szczegółowe i materiały dodatkowe oraz porady praktyczne   mogą być udzielane w ramach konsultacji

e-mail: pawula@main.amu.edu.pl

strona główna: http://main.amu.edu.pl/~pawula