Mễ HèNH TON HC TNG QUT HAI CHIU NGANG
V XM NHP MN VNG NC NGM VEN BIN
A GENERAL MATHEMATICAL MODEL
OF TWO-DIMENSIONAL HORIZONTAL FLOW
OF SEAWATER INSTRUSION IN COASTAL AQUIFERS
TRN VN MINH
UBND Thnh ph Nng
NGUYN TH HNG
Trng i hc Bỏch khoa, i hc Nng
TểM TT
Trong bi bỏo ny cỏc tỏc gi gii thiu mụ hỡnh toỏn xõy dng c tng quỏt hai chiu
ngang mụ t xõm nhp mn vo nc ngm vựng ven bin ỳng cho trng hp dũng nc
ngm cú ỏp cng nh khụng ỏp.
Thut toỏn v chng trỡnh tớnh mụ hỡnh theo phng phỏp phn t hu hn dng yu
Galerkin c thit lp d bỏo v trớ, hỡnh dng mt phõn cỏch gia nc nht v nc
mn trong cỏc iu kin khai thỏc ngun nc nht khỏc nhau.
Kt qu tớnh toỏn trờn c so sỏnh kim chng vi cỏc phng phỏp tớnh toỏn gii tớch trong
cỏc trng hp n gin.
ABSTRACT
In this paper, the authors present a general mathematical model of two dimensional horizontal
flow of seawater intrusion into coastal confined and unconfined aquifers. Algorithms and
program of this model are formulated by weak Galerkin finite element method for prediction of
the transient effect of pumping well on seawater instrusion into coastal confined and unconfined
aquifers. The validity of the model is tested by using the analytical solutions.
1. Mở đầu
Vùng ven bờ biển, nguồn nớc nhạt chảy ra biển, hình thành vùng chuyển tiếp giữa nớc
nhạt và nớc mặn. Trong vùng này mật độ của nớc hỗn hợp biến đổi từ nớc nhạt đến nớc mặn.
Tuy nhiên dới những điều kiện nhất định, chiều rộng của vùng này là nhỏ so với chiều dày của
tầng chứa nớc; vì vậy vùng chuyển tiếp từ nớc nhạt đến nớc mặn có thể xem là một mặt phân
cách. Khi nhu cầu khai thác nớc ngầm tăng lên làm giảm lợng nớc nhạt đổ về biển, mặt phân
cách sẽ di chuyển vào phía trong đất liền. Hiện tợng này đợc gọi là sự xâm nhập của nớc biển.

a
)h + .(T
a
f) = -I
s
- q'
s
(2)
Trong đó:
f = h
f
/ ; T = K(H
1
+h
f
) ; T
a
= T(H
1
- h)/(H
1
+ h
f
)
Hệ phơng trình chỉ đạo dòng chảy xâm nhập mặn của nớc ngầm vùng ven biển hai
chiều ngang có áp [3] đợc cho nh sau:
-.(T f) + .(T
a
h) = I
f

)
ở đây T đợc gọi là hệ số dẫn truyền của tầng thấm; đối với tầng thấm không đẳng hớng ta có T
= (T
x
,T
y
)
Từ các hệ phơng trình (1), (2) và (3), (4) ở trên đã thiết lập riêng lẻ cho dòng chảy xâm
nhập mặn hai chiều ngang không áp và có áp; chúng tôi viết lại hệ phơng trình nầy ở dạng
tổng quát đúng cho cả dòng chảy không áp lẫn có áp nh sau:
'
s
'
fsfayaxyx
qqII)
y
h
T(
y
)
x
h
T(
x
)
y
f
T(
y
)





+






(5)
,
ssayaxayax0
qI)
y
f
T(
y
)
x
f
T(
x
)
y
h
T(
y
)







+








(6)
ở đây:
f
fsff
;
h)1(
f


=

+
=
[ ] [ ]
2f1yy2f1xx

S
0
- độ trử nớc riêng ; t- thời gian quan sát
h
f
- chiều cao mực nớc ngầm so với mực chuẩn

f
,
s
- cột nớc thuỷ lực của vùng nớc nhạt, nớc mặn tơng ứng
I
f
, Is- lợng nớc cung cấp cho nớc nhạt, nớc mặn trong tầng thấm nớc
q'
f
, q'
s
- các điểm nguồn hoặc điểm tụ trong vùng nớc nhạt và nớc mặn
Kx, Ky- hệ số thấm theo phơng x, y
Hệ phơng trình vi phân ở trên bao gồm các phơng trình (5) và (6) là những phơng trình
vi phân đạo hàm riêng phi tuyến, việc giải tìm nghiệm chính xác của hệ phơng trình rất khó
thực hiện.
3. Phơng pháp giải số
Để giải gần đúng các hệ phơng trình trên chúng tôi sử dụng phơng pháp phần tử hữu
hạn Galerkin đợc thiết lập ở dạng chuẩn, với phần tử tam giác và giả sử rằng tại mỗi phần tử
các giá trị của hệ số thấm theo các phơng x, y là không thay đổi ta có các phơng trình thiết lập
cho một phần tử nh sau:

Trong đó:

Sau khi thực hiện các phép biến đổi toán học, ta nhận đợc hệ phơng trình tuyến tính
sau:
[ ] [ ]














+
+
+






































n
1
2f1

A
)II(
.
.
3
A
)II(
h
.
.
h
K
f
.
.
f
K
(9)
)7()( dAPLN
m
A
i

[ ] [ ]
[ ] [ ]




















=
















1
)e(
1s
tt
n
1
b
n
1
b
h
.
.
h
t
C
q
.
.
q
3
AI
.
.
3
AI
h
.
.
h

)e(
sisi
m
1e
)e(
si
)e(
fisifi
m
1e
m
1e
)e()e(
bb
m
1e
)e(
aa
===+=+
===


===
= ==
[ ]






i
ay
)e(
2
kjkik
kj
2
jij
kiji
2
i
ax
)e(
)e(
b
ccccc
ccccc
ccccc
A4
T
bbbbb
bbbbb
bbbbb
A4
T
K
Trong đó:
So
(e)
- độ trử nớc của phần tử tam giác













+













=
100
010
001

2
jij
kiji
2
i
)e(
x
)e(
a
Bờ
biển
y
x
Q
1
Q
2
q
n1
q
n2
A
B
C D
Tại thời điểm ban đầu t = 0 chiều sâu (h) của mặt phân đợc xác định theo công thức
của Ghyben-Herzberg:
f
fs
f
hh

i
đã tính đợc ở bớc thời gian trớc.
4. Kiểm tra mô hình
4.1. Số liệu tính toán
Khu vực tính toán là Hoà Khánh, thành phố Đà Nẵng, khi có ba giếng bơm khai thác
hoạt động theo sơ đồ không áp Hình 1(a) với các số liệu sau:
Hệ số thấm K= 8.1 m/ngày; lợng nớc ngấm đều I = 0 m/ngày; độ trữ nớc riêng So =
0.6; chiều sâu của tầng không thấm Ht = 25 m; lu lợng các giếng bơm tại các điểm có toạ độ
(x,y): Q
1
(x
1
,y
1
)=Q
1
(1000,0)=120m3/ngày, Q
2
(x
2
,y
2
)=Q
2
(1000,100)= 60m3/ngày, Q
3
(x
3
,y
3