1
ĐẠI HỌC QUỐC GIA HÀ NỘI
TRƯỜNG ĐẠI HỌC KHOA HỌC TỰ NHIÊN
Nguyễn Đức Hạnh
NGHIÊN CỨU XÂY DỰNG MÔ HÌNH THỦY ĐỘNG LỰC BA
CHIỀU TÍNH TOÁN TRƯỜNG DÒNG CHẢY XUNG QUANH CÁC
CÔNG TRÌNH THỦY LỰC PHỨC TẠP LUẬN VĂN THẠC SĨ KHOA HỌC
Hà Nội – Năm 2013
2 ĐẠI HỌC QUỐC GIA HÀ NỘI
TRƯỜNG ĐẠI HỌC KHOA HỌC TỰ NHIÊN
2.2.1. Mô hình nguyên gốc của Hosoda 24
2.2.2. Chỉnh sửa mô hình để có thể chạy với trường hợp kè hoàn lưu và
kè chảy ngập 30
2.2.2. Sơ đồ khối mô hình 31
2.2.3. Chương trình tính toán 32
Chương 3 – MỘT SỐ KẾT QUẢ THỬ NGHIỆM MÔ HÌNH VỚI CÁC
CÔNG TRÌNH THỦY LỰC PHỨC TẠP 43
3.1. Thử nghiệm mô hình với thí nghiệm đoạn sông thẳng có công trình 43
3.1.1. Thử nghiệm với thí nghiệm số 43
3.1.2. Thử nghiệm với thí nghiệm vật lý của Tominaga 47
3.2. Thử nghiệm mô hình với thí nghiệm vật lý đoạn sông cong có công trình
55
3.2.1. Thí nghiệm vật lý 55
4
3.2.2. Kết quả kiểm nghiệm 57
KẾT LUẬN 62
TÀI LIỆU THAM KHẢO 63
PHỤ LỤC 66
5
DANH MỤC CÁC HÌNH
Hình 1. Mô tả thí nghiệm của Munita và Shimizu (1994) trên mặt bằng (hình có
tính minh họa không chính xác về tỷ lệ) 26
Hình 2. Minh họa lưới tính toán sử dụng trong mô hình 27
Hình 3. Trường vận tốc theo mặt cắt dọc A-A thượng lưu kè mỏ hàn 28
Hình 4. Trường vận tốc trên mặt cắt ngang B-B 28
Hình 5. Trường vận tốc trên mặt cắt ngang C-C 28
Hình 6. So sánh lưu tốc tính toán và thực đo trên các mặt cắt ngang a- Mặt B-B
b- Mặt C-C c- Mặt D-D 29
Hình 25. Cấu tạo tấm hướng dòng 56
Hình 26. Mặt bằng lưới tính và vị trí các công trình. 57
Hình 27. So sánh mực nước giữa thí nghiệm và mô phỏng trong trường hợp có
công trình tại các mặt cắt 58
Hình 28. So sánh phân bố vận tốc (Uave = u2 + v2) trong trường hợp có công
trình và không có công trình 59
Hình 29. So sánh về sự phân bố vận tốc tốc (Uave = u2 + v2) giữa thí nghiệm vật
lý và mô phỏng bằng mô hình trong trường hợp có công trình 60
Hình 30. Mặt bằng các véc-tơ vận tốc ở trên mặt và dưới đáy 61
7
MỞ ĐẦU
Giao thông đường thủy trên các hệ thống sông chính ở Việt Nam đóng vai trò
rất quan trọng trong việc vận chuyển hàng hóa ở miền đồng bằng. Tuy nhiên giao
thông thủy lại phụ thuộc rất nhiều vào lưu lượng dòng chảy, độ sâu dòng chảy và
địa hình đáy sông. Việc xác định luồng chính trên các hệ thống sông lớn là rất cần
thiết và để tạo thuận lợi cho các phương tiện tham gia giao thông thủy, phòng tránh
tai nạn cũng như các thiệt hại có thể xảy ra thì cần phải duy trì sự ổn định của các
luồng chính. Để làm được điều đó, trong thực tiễn đã có một số các công trình kè
mỏ hàn, kè hướng dòng nhằm mục đích điều khiển dòng chảy tại một số vị trí và
làm thay đổi các quá trình vận chuyển trầm tích của dòng chảy sông theo hướng tạo
các luồng chính ổn định và bồi lấp các luồng phụ nhằm tăng độ sâu dòng chảy trong
mùa kiệt phục vụ giao thông thủy thuận lợi.
Mặt khác, nhằm mục đích bảo vệ bờ tại những khu vực trọng điểm xói lở
người ta thường sử dụng các loại công trình như kè lát mái, đập mỏ hàn, hoặc nhằm
mục đích bảo vệ các chân trụ cầu thì dạng công trình thường gặp đó là các mố trụ
cầu có hình dạng gần giống các đập mỏ hàn. Để mô tả các loại công trình có cấu
trúc tương tự như vậy, trong thủy lực công trình thường dùng cụm từ công trình kè
mỏ hàn và tương tự (spur-dike-like structures).
Tuy nhiên, một nguyên nhân thường gặp gây nên sự phá hoại công trình là
lực phức tạp”, tập trung vào nghiên cứu tìm hiểu một mô hình số trị ba chiều đã
được xây dựng bởi GS. Hosoda và các cộng sự [20], từ đó thay đổi lại mã nguồn
chương trình gốc để làm cho chương trình gốc có thể mô tả được với những trường
hợp có thể xuất hiện trong thực tiễn ở Việt Nam. Mô hình được kiểm nghiệm lại với
một số các kết quả thí nghiệm. Bố cục của luận văn bao gồm:
Mở đầu
Chương 1: Tổng quan các nghiên cứu mô hình ba chiều tính toán trường
dòng chảy
Chương 2: Xây dựng mô hình thủy động lực ba chiều tính toán trường dòng
chảy xung quanh các công trình thủy lực phức tạp.
Chương 3: Một số kết quả thử nghiệm mô hình
Kết luận
Tài liệu tham khảo
9
Chương 1 - TỔNG QUAN CÁC NGHIÊN CỨU MÔ HÌNH
BA CHIỀU TÍNH TOÁN TRƯỜNG DÒNG CHẢY
1.1. Tổng quan tình hình nghiên cứu trường dòng chảy ba chiều xung quanh
các công trình thủy lực
Để nghiên cứu, đánh giá được quá trình bồi xói xung quanh các công trình
thủy lực trên sông như trụ cầu, kè đập, mỏ hàn,… thì đầu tiên phải nghiên cứu, tính
toán được trường dòng chảy xung quanh các công trình thủy lực đó. Hiện nay, mô
hình toán và mô hình vật lý là hai phương pháp quan trọng trong nghiên cứu trường
thủy động lực xung quanh các công trình thủy lực.
Trong những thập kỷ gần đây, các vấn đề liên quan đến xói cục bộ xung
quanh các công trình thủy lực trên sông đã được nghiên cứu chủ yếu bằng việc xây
dựng các mô hình vật lý (Nagata và cộng sự, 2005) [43] của kè mỏ hàn, chân đế trụ
cầu, trụ cầu… Các nghiên cứu về chân đế trụ cầu và kè mỏ hàn có thể được xem là
tương tự nhau do thông thường thì chúng có hình dáng gần giống nhau.
Altinbilek 1991; Kothyari và cộng sự 1992 [27]; Melville và Chiew 1999) [38].
Một danh mục rất nhiều các nghiên cứu liệt kê trên đây đã chỉ ra tầm quan
trọng của vấn đề xói lở cục bộ sau các công trình, và tất cả đều đi đến kết luận rằng,
bằng cách tiếp cận lý thuyết hay thực nghiệm, cấu trúc dòng chảy xung quanh các
công trình là hết sức phức tạp và hình dạng hố xói cục bộ phụ thuộc rất nhiều vào
các ảnh hưởng của dòng chảy, các đặc trưng bùn cát, kết cấu nền đáy cũng như hình
dạng của công trình. Tổng quan chi tiết hơn về các nghiên cứu đó có thể được tìm
thấy trong công trình của Breuser và nnk (1977) hoặc của Dey (1997) [10].
Ngày nay, động lực học chất lỏng tính toán đã trở thành một công cụ hữu
hiệu trong phân tích và thiết kế các công trình thủy lực song hành cùng với các quan
trắc thực địa và thí nghiệm trong phòng. Dẫu rằng đã có những bước tiến dài và có
11
nhiều thành tựu đã đạt được trong lĩnh vực này, tuy nhiên, việc mô phỏng số trị
trong hầu hết các vấn đề liên quan đến công trình thủy lực trong thực tế vẫn còn đặt
ra nhiều thách thức cần vượt qua thậm chí là cả với những công cụ thủy động lực
tiên tiến nhất hiện nay. Đó là do trong hầu hết các trường hợp quan tâm thì dòng
chảy diễn ra trong các miền phức tạp hoặc là nhân tạo hoặc là tự nhiên, có cấu trúc
ba chiều rõ ràng và là dòng chảy rối thực sự. Dòng chảy đó lại được thống trị bởi
các xoáy liên kết, quy mô lớn và bất ổn định. Những khó khăn trong việc mô phỏng
số trị những dạng dòng chảy đó lại còn được chồng chất thêm do những sự khác
biệt về quy mô giữa các công trình nhân tạo quy mô nhỏ (vốn làm xuất hiện các đặc
tính thủy lực phức tạp cục bộ) với quy mô lớn của các đoạn sông tự nhiên chứa các
công trình đó (vốn xác định các đặc trưng dòng chảy đến khu vực có công trình).
Một điều kiện tiên quyết để động lực học chất lỏng tính toán đóng vai trò
quan trọng trong mô hình hóa các quá trình xói lở đó là phải có một phương pháp số
trị có khả năng dự báo chính xác bức tranh thủy động lực phức tạp, bất ổn định và
được tăng cường do sự có mặt của các mố và trụ cầu. Các phương pháp đó cần phải
(1) giải quyết được một cách hiệu quả và chính xác các chi tiết phức tạp của địa
hình các mố trụ cầu thực; (2) tính toán được các địa hình quy mô lớn của đoạn sông
không thích hợp cho các mô phỏng dòng chảy xung quanh các chướng ngại nơi mà
dòng chảy ba chiều chiếm ưu thế. Mặc dầu vậy, các mô hình rối đó đã được sử dụng
trong các nghiên cứu của Michiue và Hinokidani (1992); Jia và Wang (1993); Olsen
và Malaaen (1993); Jia và Wang (1996); Ouillon và Dartus (1997); Olsen và
Kjellesvig (1998); Richardson và Panchang (1998); Sinha và cộng sự (1998) [57];
và Neary và cộng sự (1999).
Để khắc phục các tồn tại trên cần sử dụng hệ phương trình Navie-Stokes viết
ở dạng trung bình Reynolds ba chiều, đầy đủ với giả thiết phân bố áp suất phi thủy
tĩnh được khép kín với mô hình rối không cân bằng k- để mô phỏng trường dòng
chảy (Nagata và nnk, 2005). Điều này cũng được khẳng định trong nghiên cứu của
13
Ge và Sotiropoulos (2005) [17], Ge và cộng sự (2005) [18], Olsen (2003), Lai và
cộng sự (2003a, b) [30, 31], Wu và nnk (2000), Sinha và nnk (1998),
Bên cạnh đó, nhìn chung thì địa hình xung quanh các công trình thường có
dạng rất phức tạp và để mô tả chi tiết được trong hệ tọa độ Đề các cần phải chia lưới
tính rất mịn, nhưng điều đó sẽ làm giảm khả năng ứng dụng của các mô hình vì nó
yêu cầu thời gian tính toán rất lớn do việc tăng số ô lưới đồng thời giảm bước thời
gian tính. Một cách khắc phục đơn giản đó là sử dụng hệ lưới vuông góc biến đổi, ở
các khu vực không có công trình hoặc địa hình không có những biến đổi bất thường
thì độ phân giải của lưới có thể thô hơn, và ở những khu vực xung quanh công trình
hoặc nơi địa hình có thay đổi mạnh sẽ sử dụng lưới mịn hơn. Các nghiên cứu theo
hướng này đã có nhiều ứng dụng tuy nhiên, do vẫn sử dụng lưới vuông góc thông
thường nên các nghiên cứu đó được triển khai cho các mô hình trong phòng thí
nghiệm, khó áp dụng cho các đoạn sông thực tế nhất là các đoạn sông cong. Các mô
hình sử dụng lưới cong (trực giao và không trực giao) được đề xuất để giải quyết
khó khăn trên và cho phép mô tả dòng chảy trong các đoạn sông cong, nhưng với
các mô hình vận chuyển bùn cát và bồi xói thì biên của các miền tính toán luôn thay
đổi do diễn biến lòng và bờ sông. Và để khắc phục thì Nagata và các cộng sự (2005)
đã sử dụng một hệ lưới di động để thích hợp với cả bề mặt đáy và mặt nước tự do.
động. Tuy nhiên để hiểu khả năng, định lượng tính chính xác và thiết lập thí nghiệm
cho việc mô phỏng nó cần kiểm tra những mô hình này với nhiều trường hợp
nghiên cứu khác nhau.
Những khó khăn trong việc đòi hỏi số liệu phân giải không gian tốt có chất
lượng cao trong các sông đã dẫn tới chủ yếu các xây dựng và kiểm nghiệm mô hình
số hiện mới chỉ tập trung vào nghiên cứu những địa hình sông được đơn giản hóa
trong các môi trường được kiểm soát. Việc sử dụng những mô hình phòng thí
nghiệm lược bỏ đi những sự phức tạp sự có mặt độ nhám biên ví dụ như có thực vật
cũng như tạo thuận lợi cho việc thu thập số liệu có chất lượng cao.
Những nghiên cứu số ba chiều nhìn chung đều tập trung vào các kênh lăng
15
trụ hình chữ nhật có thành nhẵn, được xây dựng sử dụng các mô hình rối bậc cao
quan hệ với mô hình rối k-ε chuẩn và không có hiệu chỉnh thành phần vận tốc theo
phương thẳng đứng. Những nghiên cứu này gồm có Sofialidis và Prinos (1998),
Cokljat và Younis (1995), Thomas và Williams (1995) và một số tác giả khác đã
mô phỏng kênh bất đối xứng thẳng có thành nhẵn của Tominaga và Nezu (1991).
Những mô phỏng này gồm có kênh thẳng có địa hình đơn giản hóa: độ dốc thành
bên thẳng đứng và đáy nằm ngang. Leschziner và Rodi (1979), Shimizu và nnk
(1990), Ye và McCorquodale (1998) đã mô phỏng dòng chảy trong các kênh cong
mặt cắt ngang hình chữ nhật hoặc hình thang. Leschziner và Rodi (1979) đã mô
phỏng một khúc cong 1800 (Rozovskii 1965) và so sánh sự dâng mặt nước dự đoán
với những kết quả đo đạc cho thấy khá phù hợp. Tuy nhiên, so sánh vận tốc điểm có
hạn chế do biểu diễn của Rozovskii về vận tốc hướng ngang. Thành phần chênh
lệch của vận tốc hướng ngang tỷ lệ với giá trị trung bình độ sâu tương ứng, đã được
lấy bằng giá trị vận tốc hướng ngang, và vì vậy đây không phải là những giá trị
tuyệt đối. Tương tự như vậy Ye và McCorquodale (1998) đã mô phỏng một đoạn
cong 270
0
(Hicks 1985) có mặt cắt ngang nửa hình thang (bờ ngoài có một cạnh dốc
cấp. Hơn nữa, nghiên cứu của Sinha và nnk (1998) không đo đạc thành phần vận
tốc theo phương thẳng đứng.
Một nghiên cứu giải số ba chiều đoạn sông uốn khúc có mặt cắt ngang mang
tính thực tế (thành bên và đáy không phẳng) mà tập trung vào so sánh giữa các vận
tốc đo đạc và tính toán đã được đưa ra bời Demuren (1993). Trường hợp dòng chảy
này được nghiên cứu thực nghiệm bởi Almquist và Holley (1985). Các cao độ đáy
được đo đạc sử dụng các thiết bị đo điểm. Cường độ vận tốc cơ bản được đo bằng
một ống Pitot tĩnh, và góc dòng chảy có một cánh quạt đồng. Trong nghiên cứu này
đã có những cải thiển đáng kể về mặt chất lượng số liệu. Demuren cũng tập trung
vào so sánh các vận tốc hướng ngang và hướng dọc trung bình độ sâu, hạn chế
trong việc tìm hiểu trường dòng chảy ba chiều.
Trong nghiên cứu của Wilson và nnk (2003) đã đưa ra sự kiểm chứng về một
mô hình số ba chiều khép kín rối k-ε chuẩn, cùng với độ phân giải không gian cao,
tập số liệu ba chiều chất lượng cao (Boxall 2000). Số liệu này được thu thập từ một
kênh uốn khúc trong phòng thí nghiệm quy mô lớn có địa hình giả tự nhiên phức
tạp. Kênh thí nghiệm được tự hình thành và uốn khúc, xây dựng trong cát giới hạn
17
trong các biên bê tông và sau đó được cố định bằng việc làm rắn hóa học (Guymer
và nnk 1999). Theo nhóm tác giả thì đây là một trong những tập dữ liệu thí nghiệm
chi tiết về mặt không gian nhất có hình dạng hợp lý về mặt vật lý, và vì vậy cung
cấp một tập số liệu tiêu chuẩn tốt cho việc kiểm tra mô hình số.
Như vậy, qua phân tích tổng quan trên đây cho thấy, cho dù đã có rất nhiều
các nghiên cứu về mô hình số trị ba chiều trong đoạn sông cong cũng như các
nghiên cứu về tác động của công trình chỉnh trị, nhưng hiện vẫn chưa có mô hình
thủy động lực 3 chiều thể hiện được tính phức tạp của dòng chảy xung quanh khu
vực công trình và có tính ứng dụng đối với các công trình dạng kè mỏ hàn ngập và
kè hoàn lưu như thực tế Việt Nam hiện nay, đặc biệt là đối với các công trình thủy
lực được đặt ở vị trí các đoạn sông cong. Để đáp ứng được nhu cầu đó, nghiên cứu
này đã cố gắng cải tiến và xây dựng bộ mô hình mô phỏng trường dòng chảy ba
J
U
j
j
(1)
Phương trình động lượng:
j
j
j
G
jij
G
j
j
i
tJJ
UU
J
UUU
J
U
t
u
u
mllm
jl
i
m
j
j
j
ij
uu
xx
JJ
p
J
g
g
j
j
i
i
u
x
U
j
G
j
i
i
G
m
j
i
m
m
i
j
m
ij
u
x
u
x
19
và:
t
: thời gian;
321
: thành phần phản biến của vận tốc lưới;
u
:
véc tơ vận tốc
321
,, uuuu
;
i
u
: thành phần véc tơ vận tốc trong hệ tọa độ Đề các;
i
G
u
: thành phần của véc tơ vận tốc lưới trong hệ tọa độ Đề các;
f
: véc tơ trọng lực
g,0,0f
;
g
: gia tốc trọng trường;
ij
: tensor ứng suất nhớt; và
ji
uu
J
U
J
U
J
U
zyx
(3)
hoặc
0
1
11
(4)
Phương trình động lượng:
phương –x
zyx
g
J
p
g
p
g
p
g
JzyxtJ
wvu
zyx
wvu
G
zxy
G
yxx
G
x
x
,,,0,0
11
,,
,,
,,,,
,,,,
,,,,20
yz
wu
xz
vw
zy
vv
yy
vu
xy
uw
zx
uv
yx
uu
xx
ww
zz
wv
yz
wu
xz
vw
zy
vv
yy
vu
xy
uw
zx
uv
zxyxxx
(5)
phương –y
J
UUU
J
U
t
yzyyyx
z
G
z
y
G
y
x
G
x
yz
G
zyy
G
yyx
G
x
y
''''''
''''''
''''''
''''''
''''''
''''''
''''''
''''''
''''''
1
ww
zz
wv
yz
wu
xz
vw
zy
vv
yy
vu
xy
uw
zx
uv
yx
vu
xy
uw
zx
uv
yx
uu
xx
J
zzyzxz
yzyyxy
zxyxxx
zzyzxz
yzyyxy
zxyxxx
zzyzxz
yzyyxy
zxyxxx
wvu
zyx
wvu
J
UU
zyy
wvu
J
UU
zyx
wvu
J
UU
J
UUU
J
UUU
J
UUU
J
U
t
zzzyzx
z
G
z
y
G
y
x
''''''
''''''
''''''
''''''
''''''
''''''
''''''
''''''
''''''
xy
uw
zx
uv
yx
uu
xx
ww
zz
wv
yz
wu
xz
vw
zy
vv
yy
vu
xy
uw
zx
uv
yx
uu
xx
J
zzyzxz
yzyyxy
zxyxxx
zzyzxz
bãi) vốn có quy mô thời gian lớn hơn là các thành phần dao động rối, một mô hình
k- phi tuyến chứa các hàm thực nghiệm đã được sử dụng và cho thấy được khả
năng mô phỏng tương đối tốt (Hosoda và cộng sự 1999, Kimura và Hosoda 2003,
Kimura và nnk 2004). Mô hình rối phi tuyến này, bao gồm thành phần bậc hai, đã
được chứng tỏ là tương đương với mô hình nhớt rối hiện phi tuyến trong mô hình
ứng suất Reynold đại số (Pope 1975, Gastki và Speziale 1993).
Năng lượng rối động học
k
và hệ số khuếch tán rối
nhận được từ các
phương trình vận chuyển sau đây trong hệ tọa độ khớp biên di động:
k
t
ll
i
j
lji
j
G
j
j
k
J
g
D
J
u
x
J
uu
J
kUU
J
k
t
Jt
2
2
''
(9)
trong đó:
00,1
k
,
30,1
,
44,1
l
C
3
1
3
1
3
2
''
ijij
t
ij
m
j
i
m
m
1
r
i
j
r
r
j
i
r
ij
x
t
và các thông số C
µ
, C
1
, C
2
, C
3
có thể được biểu diễn theo
Nagata và nnk (2005) như sau:
2
k
CD
t
2
09,01
3,0
24
;max SM
2
2
1
i
j
i
j
j
i
x
u
x
uk
2.2. Thuật toán và phương pháp giải
2.2.1. Mô hình nguyên gốc của Hosoda
a) Rời rạc hóa hệ phương trình – sơ đồ sai phân.
Phương trình động lượng và phương trình chuyển động của k và ε được rời
rạc hóa bằng phương pháp thể tích hữu hạn trên hệ thống lưới so le. Để thỏa mãn
tính liên tục cục bộ, nghiên cứu này sử dụng phương pháp Mark-and-Cell đơn giản
hóa (HSMAC) (Hirt and Cook 1972) vì cả vận tốc và áp suất có thể tính đồng thời
trong khi trong phương pháp MAC nguyên bản (Harlow and Welch, 1965) thì
phương trình Poission của biến áp suất phải được giải trước. Yếu tố đối lưu được
rời rạc hóa bằng phép nội suy ngược bậc 2 (quadratic upstream interpolation) sơ đồ
động học đối lưu QUICK (Leonard, 1979) và yếu tố khuếch tán được rời rạc hóa
bằng sơ đồ sai phân trung tâm. Trong khi xác định vận tốc cho các bước thời gian
tiếp theo khi đã biết vận tốc và áp suất cục bộ, phương trình liên tục được sử dụng
để tính áp suất thông qua phép lặp số (SOLA), đồng thời giá trị mực nước cũng
được xác định theo điều kiện động lực ở mặt thoáng. Khi giá trị mực nước này được
xác định, thì sơ đồ lưới trong mặt phẳng thẳng đứng được xây dựng bằng cách biến
đổi khoảng cách lưới. Lưu ý không có sự thay đổi lưới trong mặt phẳng (x
1
Rodi 1982; Gibson và Rodi 1989; Naot và cộng sự. 1993) được đưa ra để tính tác
động tắt dần của hiện tượng chảy rối do sự xuất hiện của các mặt thoáng. Vận tốc
trượt xác định từ profin vận tốc loga được sử dụng để xác định k
b
và ε
b
. Trong
nghiên cứu hiện nay, một số phương pháp tính toán được Kimura và Hosoda sử
dụng: (1) độ nhớt rối được biến đổi bằng hàm tắt dần Isawa và Hosoda đưa ra
(1990); và (2) tốc độ tiêu tán ở mặt thoáng được xác định bằng công thức Siguyama
và cộng sự đề xuất (1997), tương tự với quan niệm của Naot và Rodi (1982). Rõ
ràng cường độ rối theo hướng thẳng đứng giảm theo hướng mặt thoáng. Để tính
toán tác động này, nhân độ nhớt rối với hàm tắt dần f
d
được Iwasa và Hosoda đưa ra
(1990):
= 1
3
3/2